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Preface 


Nuclear thermal hydraulics is the application of thermofluid mechanics within the nuclear indus- 
try. Thermal hydraulic analysis is an important tool in addressing the global challenge to reduce 
the cost of advanced nuclear technologies. An improved predictive capability and understanding 
supports the development, optimisation and safety substantiation of nuclear power plants. 


This document is part of Nuclear Heat Transfer and Passive Cooling: Technical Volumes and Case 
Studies, a set of six technical volumes and four case studies providing information and guidance 
on aspects of nuclear thermal hydraulic analysis. This document set has been delivered by Frazer- 
Nash Consultancy, with support from a number of academic and industrial partners, as part of 
the UK Government Nuclear Innovation Programme: Advanced Reactor Design, funded by the 
Department for Business, Energy and Industrial Strategy (BEIS). 


Each technical volume outlines the technical challenges, latest analysis methods and future direc- 
tion for a specific area of nuclear thermal hydraulics. The case studies illustrate the use of a subset 
of these methods in representative nuclear industry examples. The document set is designed for 
technical users with some prior knowledge of thermofluid mechanics, who wish to know more about 
nuclear thermal hydraulics. 


The work promotes a consistent methodology for thermal hydraulic analysis of single-phase heat 
transfer and passive cooling, to inform the link between academic research and end-user needs, 
and to provide a high-quality, peer-reviewed document set suitable for use across the nuclear 
industry. 


The document set is not intended to be exhaustive or provide a set of standard engineering ‘guide- 
lines’ and it is strongly recommended that nuclear thermal hydraulic analyses are undertaken by 
Suitably Qualified and Experienced Personnel. 


The first edition of this document set has been authored by Frazer-Nash Consultancy, with the 
support of the individuals and organisations noted in each. Please acknowledge these documents 
in any work where they are used: 


Frazer-Nash Consultancy (2021) Nuclear Heat Transfer and Passive Cooling, 
Study A: Liquid Metal CFD Modelling of the TALL-3D Test Facility. 
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Introduction 


As part of the development of a new Nuclear Power Plant (NPP), it is necessary to ensure that any 
analysis tools used to support the design or safety justification are appropriate and adequately vali- 
dated. The level of validation required to support a nuclear safety case should be determined using 
a graded approach, and should ideally follow a hierarchical process that separates and simplifies 
the physical phenomena involved in the system of interest (Volume 1, Section 4.3.2). 


This case study considers the application of a simple hierarchical validation exercise (basic test 
and Separate Effect Test (SET)) to assess the suitability of Computational Fluid Dynamics (CFD) 
modelling for forced and mixed convection in liquid metals with and without Conjugate Heat Transfer 
(CHT). The specific objectives of this case study are to demonstrate: 


* The steps required as part of a validation exercise for a CFD application. 


How the low Prandtl number for liquid metal affects the turoulence modelling approach and 
heat transfer predictions. 


The use of higher fidelity Large Eddy Simulation (LES) to select the most appropriate turbu- 
lence modelling approach for forced and mixed convection flows. 


The application of CHT modelling to predict detailed solid component temperatures. 


The application of a sensitivity analysis to understand the uncertainty in the results due to 
uncertainty in the inlet and boundary conditions. 


The validation of thermal hydraulic analysis methods against experimental measurements, 
and assess the level of modelling necessary. 


As for all the case studies in this series, this analysis provides a ‘worked example’ of a specific 
modelling task to illustrate the modelling approaches described in the technical volumes. Therefore, 
it is recommended that this case study is read in conjunction with the following technical volumes: 


Volume 2: Convection, Radiation and Conjugate Heat Transfer 
Volume 3: Natural Convection and Passive Cooling 

Volume 4: Confidence and Uncertainty 

Volume 5: Liquid Metal Thermal Hydraulics 


Test data from the TALL-3D' facility in the Royal Institute of Technology (KTH), Stockholm and 
higher fidelity LES predictions are used to explore and exemplify how to validate Reynolds-Averaged 
Navier-Stokes (RANS) models, inform decision-making and improve confidence in analysis results. 


The name comes from Thermal-hydraulic ADS Lead-bismuth Loop with 3D flow test section (TALL-3D); operated with 
relevance to Accelerator-Driven Systems (ADSs). 
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Introduction 


Case Study Description 


This case study demonstrates the process of modelling liquid metal heat transfer from a position 
of no previous knowledge, and validating CFD models of mixed and forced convection flows as 
part of a hierarchical validation process (Figure 1.1). The simulations have been performed using 
the typical functionality available in commercial CFD codes in order to highlight their utility but also 
some of their limitations. 


Increasing Complexity of Phenomena 


Forced Convection Mixed Convection 


Conjugate Heat Transfer 
RANS and LES CFD 


Separate Effect Test (SET) 


Forced Convection 


Sensitivity Analysis 


RANS and LES CFD 


Basic Test 


Figure 1.1: Summary of techniques and test problems used in this case study. 


Heated Channel 


Modelling flow between two heated plates (infinite channel) provides a simple basic test case 
for CFD modelling of liquid metal. This test case describes how to apply customised turbulence 
modelling parameters, and allows their correct implementation to be checked against higher fidelity 
LES and existing Direct Numerical Simulation (DNS) results. This illustrates a typical process of 
establishing capability and confidence to help readers unfamiliar with CFD modelling of liquid metal. 


TALL-3D Test Section 


The understanding and confidence gained from modelling the heated channel can then be applied 
to the more complex geometry and heat transfer phenomena found in the 3D test section that forms 
part of the KTH TALL-3D liquid metal test loop, before considering modelling a reactor problem. 


The test facility was designed to provide data for validation of coupled CFD and 1D system codes. 
The main focus of the experimental investigations was the mutual feedback between natural cir- 
culation in the three-legged loop and stratification/mixing phenomena in the 3D test section, espe- 
cially the transition over time from forced to natural circulation of the loop when the pump stops, 
but heating is maintained. However, since this case study is focused only on the application of CFD 
analysis to the 3D test section, two statistically steady flow conditions have been selected: 


1. Initial forced circulation condition (pump running), where the heat transfer is dominated by 
forced convection. 
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2. Final natural circulation condition (pump off) with buoyancy driven natural circulation around 
the loop. This decreased flow rate, combined with local stratification and buoyancy driven 
flow, results in mixed convection heat transfer within the 3D test section. 


Two different models have been used in the TALL-3D validation exercise. 


Fluid only: Test section walls are not included in order to simplify the model and provide well- 
defined boundary conditions. A hybrid LES model has been solved to provide benchmark 
data that allows a clear comparison and assessment of RANS modelling methods. 


Fluid and solid model: Test section geometry with the walls included using a CHT approach. 
This enables the CFD predictions to be validated against the TALL-3D measurements, and 
assess the sensitivity of the results to more parameters. 


The following modelling approaches have been used to simulate the forced and mixed convection 
conditions within the 3D test section: 


* Steady RANS modelling using a 2D axisymmetric geometry. 
* Steady RANS modelling using a 3D geometry. 

¢ Unsteady RANS modelling using a 3D geometry. 

¢ LES modelling of a 3D geometry. 


This sequence of modelling approaches has an increasing computational cost, but also an increas- 
ing level of fidelity and ability to represent physical phenomena. The intention is to demonstrate 
what level of modelling is necessary to predict the experimental results and what is gained by each 
level of increased fidelity (or conversely, what is lost with each level of approximation). 
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2.1 


Heated Channel Flow 


Volume 5 (Liquid Metal Thermal Hydraulics) describes the main issue associated with modelling 
turbulent heat transfer in liquid metal — the low Prandtl number, Pr, of the fluid (0.001 < Pr < 
0.1). This means that the thermal conduction effects in the fluid are much more significant than in 
water or gases, and that the thermal boundary layer at a wall is much thicker than the momentum 
boundary layer. This means that the assumption that turbulent Prandtl number (Pr,) should be = 1, 
based on the Reynolds analogy, is no longer valid. 


Most RANS turbulence models use the turbulent Prandtl number concept to relate turbulent heat 
flux to turbulent momentum transfer, so it is important to understand and evaluate the impact of 
this approximation within the turbulence model on the temperature predictions for any liquid metal 
application. A simple case of fully developed turbulent flow through a uniformly heated channel has 
been used here to understand the impact of this, and review some of the correlations that have 
been proposed to replace the standard values of turbulent Prandtl number that are used (typically 
around Pr, = 0.85). This case has been selected because existing DNS and LES predictions of 
it were reviewed in the Simulations and Experiments for the Safety Assessment of MEtal cooled 
reactors (SESAME) project (Roelofs, 2019). 


Problem Definition 


Fully developed turbulent flow of a low Prandtl number fluid through a uniformly heated channel is 
simulated using periodic boundary conditions in the streamwise (x) and spanwise (z) directions, 
see Figure 2.1. 


Figure 2.1: Fully developed turbulent flow through heated channel. 
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2.2 


Heated Channel Flow S, 


The flow is characterised by the Reynolds number based on the friction velocity (Re, = pu,6/,), 
where the friction velocity 7, = \/Tv,/p, 6 is half the channel width, 7, is the time-averaged shear 
stress at the wall, uw is the viscosity and p is the density. 


A uniform and constant heat flux qy (W/m?) is applied to the top and bottom surface, which 
generates a temperature rise in the fluid of AT = 27q,,/(pcpu), where wu is the streamwise time 
and space-averaged velocity across the channel (i.e. mass flow rate W = pL,L,@) and c, is the 
specific heat capacity. 


DNS and LES calculations have been performed at Re; = 180 and 590 for Pr = 0.01 (Bricteux 
et al., 2012, Tiselj and Cizelj, 2012), and LES calculations have been undertaken at Re, = 2,000 
for Pr =0.01 and 0.025 (Duponcheel et a/., 2014). These conditions are equivalent to bulk Reynolds 
numbers (Re = pu26/) of 5,660, 22,000 and 84,000 respectively’. 


Constant fluid properties have been used in the reference DNS and LES calculations, so that only 
forced convection is simulated and there are no buoyancy effects associated with the flow. In reality, 
properties may change with temperature and at low velocities natural convection near the wall may 
become significant. 


The results of importance from the simulation are the time-averaged velocity and temperature 
profiles across the channel, which are normalized using: 


eas 
Uy 
Tw — T he w 
gt = (Tw =T) where the friction temperature 6, = a 
0, Peps 
4 purd 
—— where d = y for y < 6 and d = (26 — y) fory > 6 


Planning the Analysis 
The purpose of the heated channel study is to: 


* Model simple liquid metal flows using a commercial CFD code. 

* Validate RANS turbulence models against reference solutions. 

« Understand the impact of Prandtl number on the temperature profile. 

* Extend the Prandtl number range up to 0.035 (highest value within TALL-3D test case). 
* Investigate the sensitivity of the temperature profile to the turbulent Prandtl number. 


The ANSYS Fluent commercial CFD software has been selected for the heated channel study, and 
subsequent TALL-3D analysis (Section 4), as there is: 


* A high level of expertise and experience in this code within the team. 
* Broad range of RANS turbulence models and LES methods available. 
* Ability to modify and customise the turbulent Prandtl number. 


This definition of bulk Reynolds number has been used as it is consistent with the previous work. However, the Reynolds 
number for an infinite channel is normally defined as Re = puDp/p, where Dp = 4A/P = 46. 
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2.2.1 


Pete 


Heated Channel Flow 


Case Selection 


In order to gain confidence in the RANS turbulence model predictions, it is necessary to validate 
the results against experimental measurements or benchmark simulations. Therefore, the cases 
in the heated channel study have been selected to match the existing available benchmark data 
(DNS and LES results). This enables the RANS turbulence models to be validated over a range of 
Reynolds numbers (Re, = 180, 590, 2,000) and liquid metal Prandtl numbers (Pr = 0.01, 0.025). 


The Prandtl number range has also been extended above 0.025 as part of this work to enable a 
more direct comparison to the TALL-3D test conditions. This requires additional LES results to be 
generated in order to validate and assess the RANS simulations at Pr = 0.035. 


In order to simulate these cases in a CFD model, the geometry, material properties, flow and 
thermal parameters need to be defined. The material properties were set to representative values 
for liquid metals to achieve the required Prandtl number (Table 2.1), as used in the reference DNS 
and LES simulations. These do not represent any real fluid, and have been chosen to give the 
appropriate Prandtl number. 


A value of 6 = 0.01 m has been used to generate the geometry with 7 set to 3.2 for simplicity, and 
the heat flux (qy) has been calculated to give AT = 2K. These values do not affect the results 
of importance as the simulation outputs are non-dimensionalised, and so representative physical 
values have been used. The list of cases selected and input parameters based on these values is 
summarised in Table 2.2. 


Prandtl number 0.01 0.025 0.035 


Density (kg/m?) 10000 10000 10000 
Specific heat capacity (J/kg K) 100 100 140 
Thermal conductivity (W/m kK) 10 10 10 
Viscosity (kg/m s) 0.001 0.0025 0.0025 


Table 2.1: Thermophysical material properties for the heated channel study. 


Case 1 2 3 4 5 6 7 


Friction Reynolds number 180 590 2,000 2,000 180 590 2,000 
Prandtl number 0.01 0.01 0.01 0.025 0.035 0.035 0.035 


Mass flow (kg/s) 0.18112 0.704 2.688 6.72 0.4528 1.76 6.72 


qu (W/m?) 8843.75 34375 131250 328125 30953.125 120312.5 459375 


Existing data DNSLES DNSLES LES LES N/A N/A N/A 


Table 2.2: Case selection for the heated channel study. 


Sensitivity to Turbulent Prandtl Number 


Specific RANS models for turbulent heat transfer in low Pr flows are the subject of current research 
and development. Advanced models, such as the Algebraic Heat Flux Model (AHFM), which are 
discussed in more detail in Volume 3 (Natural Convection and Passive Cooling), Volume 5 and 
by Roelofs (2019), are not commonly implemented and used in most commercial CFD codes. 
Therefore, this study is limited to the more commonly used turbulent Prandtl number approach and 
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2.3 


2:01 


Heated Channel Flow 


investigates the sensitivity of the temperature predictions to three proposed modifications to the 
turbulent Prandtl number in a conventional two-equation Eddy Viscosity Model (EVM). 


1. Reynolds modification: A single value of turbulent Prandtl number was proposed by Reynolds 
(1975) as a function of the bulk Reynolds number and Péclet number of the flow. 


1 


Pre = (1+ 100Pe~1/?) (, ae 
(= 


_ 0.15) where Pe = PrRe 


2. Kays modification: Kays (1994) proposed a correlation that calculates local values of Pr; as 
a function of the turbulent Péclet number, Pe; = Pr(v;/v), where v; is turbulent viscosity 
and v is the kinematic viscosity. 

0.7 


ae 


3. Weigand modification: Weigand et al. (1997) proposed an alternative correlation that de- 
pends on the bulk Reynolds number. 


1 1 1 1 
= t CP, CPe,)* |1 
Pe Pag ee | a ( care) 
100 
where C = 0.3 and Prioo = 0.85 + Pr Reo 58 


These three turbulent Prandtl number modifications are identified in Roelofs (2019) and represent 
different levels of complexity from a single global value to a local value that depends on the cell 
turbulent Péclet number and bulk Reynolds number. It is worth noting that since the Reynolds and 
Weigand modifications both depend on a bulk Reynolds number of the flow, it will be difficult to 
define or calculate the modified turbulent Prandtl number for complex flows and geometries using 
these methods. 


Performing the Analysis 


The geometry and boundary conditions for the heated channel case are simple and well defined, 
and the main objective is to investigate the sensitivity to Prandtl number and RANS turbulence 
model parameters. 


Meshing Approach 


A simple hexahedral mesh was generated for the channel using the parameters detailed in Ta- 
ble 2.3, which was developed based on the recommendations in Volume 2 (Section 3.4.2.5). 


Since a ‘low-Re’ method is being used to resolve the viscous sublayer, the First Cell Height (FCH) 
was calculated using the correlation y = 2v/u, (the factor of 2 is needed because Fluent is a 
cell-centred solver) in order to achieve a y* of 1. Due to the low Prandtl number of liquid metal, the 
thermal boundary layer is thicker than the momentum boundary layer and therefore it will be over- 
resolved. In this study, the friction velocity (G, = Re,u/p65), and hence FCH, can be calculated as 
Re, is specified for each case. 


The near-wall inflation layer consists of 22 layers of cells with a linear growth rate of 1.15. The 
global mesh size was set based on a first cell aspect ratio of 20-30 with no jump in mesh size 
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Heated Channel Flow 


at the edge of the boundary layer mesh and unit aspect ratio in the main flow (Figure 2.2). For 
simplicity, the Re; = 590 mesh has been used for Re, = 180. 


Case 1,2,5,6 3,4,7 
First cell height (mm) 0.03 0.01 
Maximum Aspect Ratio 20 30 

2D mesh (Ly x L, cells) 64x 107 130x324 
2D mesh (cells) 6,848 42,120 
3D mesh (cells) 373,248 6,823,440 


Wall refined 3D mesh (cells) 8,200,192 


Table 2.3: Mesh statistics for heated channel cases. 


Figure 2.2: Heated channel mesh: Boundary layer mesh (left), refined wall mesh (right). 


In addition, a separate refined wall mesh has been developed with an aspect ratio of 2 at the wall in 
order to properly resolve the velocity fluctuations near the wall. This was generated by creating an 
initial mesh with a coarse boundary layer. 3-4 layers of cells were then selected and isotropically 
refined (i.e. each cell split into eight). This process was undertaken three times until the FCH was 
the same as the boundary layer mesh. 


Modelling Approach 


A constant heat flux was applied to the top and bottom walls and periodic boundaries at the sides. 
The inlet/outlet surfaces are periodic with a prescribed mass flow rate and upstream bulk tem- 
perature of 300K. This approach applies a body force to the fluid to drive the periodic flow to a 
target flow rate, and subtracts a mean temperature difference from the outlet temperature field to 
recycle it to the inlet. The model has been solved using the commercial ANSYS Fluent CFD soft- 
ware with the SIMPLE pressure-velocity coupling scheme, second-order discretization and default 
under-relaxation parameters. 


The realizable k-¢ with enhanced wall treatment and k-w Shear Stress Transport (SST) RANS 
turbulence models were solved without and with modifications to the turbulent Prandtl number (Sec- 
tion 2.2.2), which were implemented using a User Defined Function (UDF). Both of these models 
are ‘low-Re’ models and the y* on the walls was checked to ensure that it was below 1, and so 
suitable to resolve the viscous sublayer using such a ‘low-Re’ method (Volume 1, Section 4.5.3.3). 


LES calculations using a number of different sub-grid scale models were also run for Case 2, which 
are discussed in Section 2.3.3. The time step was specified to give a maximum Courant-Friedrichs 
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2.3.3 


Heated Channel Flow 


Lewy (CFL) number of 0.5 within the domain (Volume 3, Section 3.2.5). The instantaneous velocity 
was initialised from a RANS solution and then solved for at least 25 flow passes through the 
domain. The time-varying results were averaged over a further 70,000 time steps, and a single 
time-averaged velocity and temperature profile was generated by averaging the time-averaged 
results in the spanwise (z) direction at the streamwise midpoint (x = L,./2). 


The suitability of the LES mesh was checked by assessing the integral length scale in the initial 
RANS simulation (Volume 3, Section 3.4.2), which confirmed that the mesh was able to resolve 
over 80 % of the turbulent kinetic energy. No formal RANS mesh sensitivity has been undertaken 
as the LES mesh has been used for all RANS simulations. 


The RANS and LES model convergence was checked by monitoring the maximum velocity and 
maximum and minimum temperature on the YZ plane at the streamwise midpoint. This confirmed 
that the RANS solution had been washed through sufficiently after 25 flow passes, and that the 
time-averaged velocity and temperatures had converged. 


The 2D steady RANS solutions converged in under 10 minutes on 8 processors, while the 3D 
steady RANS and LES solutions took less than 2 hours and almost 4 days on 80 processors 
respectively for the Re, = 2,000 case. This demonstrates the large increase in solution time asso- 
ciated with running LES models. 


Discussion of Results 


Case 2 (Pr = 0.01, Re, = 590) was solved initially for the standard Dynamic Smagorinsky LES 
model for both the boundary layer and refined wall mesh in order to investigate the impact of the 
cell aspect ratio near the wall on the results. A number of hybrid LES methods available in Fluent 
have also been investigated to understand and assess their ability to resolve liquid metal heat 
transfer. These are discussed in more detail in Volume 1 (Section 4.5.3.4), and use a RANS wall 
treatment to blend between a RANS model near the wall to LES in the rest of the domain. 


* Algebraic Wall Modeled Large Eddy Simulation (WMLES). 
* Improved Delayed Detached Eddy Simulation (IDDES) SST. 
* Stress-Blended Eddy Simulation (SBES) k-w. 


The instantaneous and time-averaged contours of velocity and temperature for the WMLES solu- 
tion are shown in Figure 2.3. The instantaneous WMLES predictions show that the fluctuations in 
the temperature field are reduced due to the high thermal conductivity, compared to the fluctuations 
in the velocity field. 


Although there is some spanwise variation in the time-averaged results, the spanwise averaged 
velocity and temperature profile was checked for convergence during the transient solution. The 
time-averaged contours also demonstrate that the momentum boundary layer is much thinner than 
the thermal boundary layer, as expected for low Prandtl number fluids. 


The LES, hybrid LES, realizable k-¢ and k-w SST results have been plotted in Figure 2.4 using 
the non-dimensional parameters described in Section 2.1. The flow through the channel is forced 
(pressure driven) with no buoyancy effects, so the velocity profile is not influenced by the tempera- 
ture profile. 
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Heated Channel Flow 
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Figure 2.3: Contours on the spanwise plane at the streamwise midpoint for Case 2. 
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Figure 2.4: Non-dimensional velocity (left) and temperature (right) profiles for Case 2. 


Heated Channel Flow 


Velocity profile: The Dynamic Smagorinsky model with refined wall mesh provides the best match 
to the DNS data (Tiselj and Cizelj, 2012) as the cells near the wall have a much lower aspect 
ratio, and so the turbulent eddies in the boundary layer are properly resolved. The three hy- 
brid LES model variants resolve the turbulence in the freestream and apply different near-wall 
modelling approaches, of these the IDDES SST is slightly better, but the level of agreement 
is broadly similar. As expected, the RANS turbulence models fit the DNS data very well, as 
the viscous sublayer is resolved and the turbulence models have been developed based on 
standard fully developed flows. 


Temperature profile: The corresponding LES and hybrid LES temperature profiles all show a very 
similar trend with reasonable match to the DNS data. Although there is not much difference, 
the Dynamic Smagorinsky and WMLES models provide the best fit in this case. Since the 
near-wall mesh requirements for LES solutions are often impractical, the WMLES model 
has been selected and used to generate hybrid LES results for the remaining cases. The 
maximum 6* values for the RANS turbulence models are 20% lower than the LES results, 
which highlights the expected issues associated with modelling low Prandtl number fluids 
using a constant turbulent Prandtl number of 0.85. 


The remaining cases in Table 2.2 have been solved using the WMLES model, realizable k-¢ and 
k-w SST RANS turbulence models with the Reynolds, Kays and Weigand modification to the 
turbulent Prandtl number (Section 2.2.2). 


The results (Figure 2.5) show that the standard constant turbulent Prandtl number always leads 
to under-prediction of the maximum 6* value, while the Weigand model is consistently over- 
estimating it. The increased turbulent Prandtl number associated with the Reynolds model gen- 
erally provides a good match to the data, but depends on the Reynolds number of the flow, which 
is not known for complex three-dimensional flows. The Kays model also compares well with the 
LES predictions and only depends on the local turbulent viscosity and material properties, which 
are calculated by default in the RANS solutions. The Kays modification can therefore be more read- 
ily applied to complex flow configurations, but may not work as well in natural or mixed convection 
regimes (Roelofs, 2019). 


The k-w SST turbulence model performs better at low Reynolds numbers, while the realizable k - 
é matches the LES data better at high Reynolds numbers. This is linked to the turbulent viscosity 
ratio profile, as lower turbulent viscosity ratios generate higher turbulent Prandtl numbers in the 
Kays model, which increases the maximum @* value. 
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Heated Channel Flow 
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Figure 2.5: Non-dimensional temperature profiles for heated channel study. 
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2.4 


Heated Channel Flow 


Application of the Results 


The heated channel study has highlighted the sensitivity of the temperature profile to the turbulent 
Prandtl number in RANS turbulence models, and the challenges of modelling low Prandtl number 
liquid metal flows, even in this simple case. 


This study has demonstrated that additional care is needed to model liquid metal flows using 
commercial CFD software, and the following steps should be taken: 


« Investigate the sensitivity of the temperature predictions to the turbulent Prandtl number and 
turbulent heat flux models in RANS turbulence models. 
* Assess the sensitivity of the results of importance to different RANS turbulence models. 


+ If possible, undertake LES modelling of the relevant phenomena to provide a benchmark so- 
lution to compare and select the most appropriate RANS modelling approach for a particular 
application. 
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TALL-3D Facility Description 


TALL-3D is an experimental test facility operated by KTH, Stockholm, that was part of the European 
collaborative Research and Development (R&D) projects funded by Euratom: Thermal Hydraulics 
of Innovative Nuclear Systems (THINS) and SESAME (Kudinov and Grishchenko, 2019). It was 
designed to provide experimental data for validation of standalone and coupled thermal hydraulic 
system codes and CFD models for liquid metals. 


The TALL-3D facility consists of the primary loop, secondary cooling loop, differential pressure 
measurement system and pressurized service loop. The fluid in the primary loop is Lead-Bismuth 
Eutectic (LBE) with a maximum temperature of 460°C and pressures up to 0.7 MPa. The total 
height of the facility is about 7m and total electric power is 80 kW, including a 27 kW main heater 
and a 15 kW heater in the 3D test section. 


The primary loop (Figure 3.1) consists of three 5.83 m long vertical legs connected by two 0.74m 
wide horizontal sections (total width 1.48 m) with 27.86 mm nominal inner pipe diameter. 


Heat 
exchanger 


Electromagnetic 


OC) pump 


Main 
heater 
3D Test 


section 


Figure 3.1: Schematic of the TALL-3D primary loop. 


Flow direction in the primary loop is generally clockwise and consists of the following legs: 


* The main heater leg (left) that includes a 870 mm long pin-type electric heater in the lower 
part with expansion tank at the top. 
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3.1 


TALL-3D Facility Description ‘y, 


« The 3D leg (middle) that includes a pool-type 3D test section, which is designed to produce 
more complex 3D flow phenomena. This includes transition from forced to mixed convection 
inside the 3D test section depending on the loop characteristics. 


¢ The heat exchanger (cold) leg (right) has a 1m long counter-current double-pipe heat ex- 
changer at the top with an electric permanent magnet pump providing forced circulation up 
to 5kg/s. 


The TALL-3D facility has been used to experimentally investigate the transition from forced to 
natural circulation or vice versa over time for liquid metal at different mass flow rates and heater 
powers. This facility has been used to validate system code models of the loop behaviour, CFD 
models of the 3D test section and the coupling between them in the THINS and SESAME projects. 


3D Test Section 


The 3D test section boundary extends from the inlet flange (1,211 mm elevation) to the outlet flange 
(2,111 mm elevation) with the top of the base plate located at an elevation of 1,461 mm. The data 
monitored at the boundaries include LBE flow temperature, differential pressure and mass flow 
rate. 


Figure 3.2 shows the geometry of the 3D test section (Grishchenko ef a/., 2017), which consists of 
a 300 mm diameter by 200 mm long cylindrical vessel with a 200 mm diameter by 10 mm circular 
inner plate located 170 mm above the base. The inner plate is attached to the top of the test section 
by four separators, which are fin-shaped to minimise their impact on the flow. The flow enters at the 
bottom through a 17mm diameter pipe and exits from the top through a 27.8mm diameter pipe. 
Two 7.5 kW rope heaters are rolled jointly around the circumference of the upper part of the test 
section with a nominal length of 116mm. 


Thermal Insulation 


e Thermocouple locations ISOVER TapeLock 7300 


Outlet pipe 
Inner pipes for in-pool TCs 


Heating element 


Heating element 
clamp 


2x7.5 kW 


Lateral wall 
(ILW & OLW) 


Fin-shaped separators} 


| 
H 


Base plate 
(BP) 


Thermal Insulation 


Inlet Pipe Nano T Ultra - Calcium Silicate 


Figure 3.2: CAD geometry and layout of the 3D test section (Kudinov and Grishchenko, 2019). 


All of the 3D test section components are manufactured from stainless steel, and surrounded 
by thermal insulation to minimise external heat loss. Two types of thermal insulation (ISOVER 
TapeLock 7300 and Nano T Ultra) are used for different regions around the 3D test section and 
pipework. 


15 of 61 


TALL-3D Facility Description e, 


In total there are 159 K-type thermocouples (1 mm diameter) around the 3D test section, including 
114 inside the LBE pool, 4 on the outer surface of the insulation and 1 on the bottom of the base 
plate (Grishchenko et a/., 2017). These thermocouples (green dots in Figure 3.2) are attached to 
the following surfaces at 4 equally spaced circumferential locations: 


* Top of Base Plate (BP) at 25mm radial intervals (16-off). 

¢ Bottom of Circular Inner Plate (CIP) at 10 mm radial intervals (41-off). 

* On the Inner Lateral Wall (ILW) at 40 mm axial intervals (16-off). 

* Inside the Inner Pipe Thermocouple (IPT) tube at 25 mm axial intervals (36-off). 
¢ On the Outer Lateral Wall (OLW) at 20 mm axial intervals (40-off). 


It is worth noting that the CIP thermocouples are aligned with the outlet fins, while the BP, ILW and 
OLW thermocouples are aligned with the IPT tube, which is located midway between the outlet 
fins. 


During normal operation, cold flow enters at the base of the 3D test section and mixes because it 
has to flow radially to pass around the inner plate. The flow is then heated through the walls by the 
rope heater before exiting through the top of the 3D test section. The mass flow rate through the 
3D test section (middle leg) is driven by the circulation (which could be forced or natural) around 
the primary loop. 


As detailed in Grishchenko e¢ al. (2017), when the mass flow rate entering the test section is 
greater than 0.7 kg/s, the pool is fully mixed with the flow dominated by the strong vertical jet hitting 
the inner plate (forced convection). As the flow rate drops below 0.3kg/s buoyancy forces become 
more important and the pool starts to stratify causing the cold plume to slow down and fall as an 
annular plume around the inner rising jet (mixed convection). 
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Planning the Analysis 


The objective of the TALL-3D case study is to validate commercial CFD codes for forced and 
natural convection in liquid metal, and understand the benefits and limitations associated with 
different approaches and turbulence models. 


This case study would generally be used as part of a hierarchical validation approach (Volume 1, 
Section 4.3.2) where the heated channel study is a basic test and the TALL-3D study is a SET. 
This provides the evidence to support and justify the use of CFD analysis for a specific reactor 
application and geometry that involves similar thermal hydraulic phenomena. 


This is intended to build on the heated channel study (Section 2), and assess the ability to model 
the following key thermal hydraulic phenomena in liquid metal. These phenomena are described in 
more detail in Volume 3 with further information that is particularly relevant to liquid metal discussed 
in Volume 5. 


Forced convection: Free circular jet expansion and impingement on a perpendicular plate with 
turbulent mixing in a plenum, coupled with forced convection from a vertical heated plate. 
This is relevant to forced circulation conditions in a liquid metal reactor with heat transfer 
(e.g. the upper and lower plenums of a pool-type reactor or the downcomer flow of a loop- 
type reactor with operating recirculation pump). 


Mixed convection: Overturning negatively buoyant circular jet (e.g. a low momentum cold jet into 
a hot stratified volume) with turbulent mixing in a plenum, coupled with buoyancy driven 
convection on a vertical heated plate and thermal stratification. This is relevant to natural 
circulation during loss of power, particularly in the lower plenum of a pool-type reactor and 
also the temperature distribution in the downcomer without forced circulation. 


Heat transfer occurs by conduction and convection through the LBE fluid, conduction through the 
stainless steel vessel and surrounding thermal insulation and convection and radiation to the ex- 
ternal ambient air. 


As described in Section 3.1, the TALL-3D test section provides a good quality liquid metal valida- 
tion test case due to the large number of thermocouples present within the vessel. This enables 
the CFD predicted temperatures to be compared against measured data in order to assess the 
accuracy of the simulation and sensitivity of the results to input parameters. In addition, this facility 
has been studied using a variety of system codes, CFD codes and coupled approaches: 


¢ Grishchenko et al. (2015) details the results of the first commissioning tests for TALL-3D 
under the THINS project, together with the pre-design and pre-test analyses results, which 
includes 2D axisymmetric CFD of the 3D test section using Star-CCM+. 
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4.1 


Planning the Analysis 


The system code and CFD analyses that were undertaken within the SESAME project are 
reported in Kudinov and Grishchenko (2019). This includes a RELAP5 model of the loop 
coupled with a 2D axisymmetric Star-CCM+ model of the 3D test section. The CFD analy- 
sis included solution verification and a detailed sensitivity analysis with System Response 
Quantities (SRQs) for flow rates of 0.3kg/s and 2.0 kg/s (Jeltsov et al., 2019). 


Moreau et al. (2019) provides a summary of the lessons learnt from the SESAME project, 
which identifies an inconsistency in the global heat balance in the TALL-3D test section of 
the order of 10%. 


Papukchiev et al. (2019) details the validation of the ATHLET thermal hydraulic system code 
for the simulation of transient flows in the TALL-3D loop. 


Grishchenko et al. (2020) presents the results of an open and blind benchmark on natural cir- 
culation instability using a range of different coupled system and CFD codes. These coupled 
methods included 2D axisymmetric CFD models in TRUST, TrioCFD and ANSYS CFX. 


Modelling Approach 


The transition from forced to natural circulation generally takes over an hour to fully develop and sta- 
bilise in the TALL-3D facility; this is too computationally expensive to resolve with a CFD analysis. 
Therefore, rather than model the transition over time, it is more appropriate to start by undertaking 
a separate steady-state CFD simulation of the initial (forced convection) and final (mixed convec- 
tion) conditions in the 3D test section when the flow and temperatures are stable i.e. statistically 
steady flow. 


Most previous CFD analyses have employed a 2D axisymmetric approach because the geometry 
of the 3D test section is circumferentially uniform, except for the presence of the four IPT tubes and 
fin-shaped separators, bolts and thermocouple mounting tubes. However, in this case study it was 
decided to undertake the validation from first principles using a sequential approach, to determine 
what level of modelling approximation is appropriate. 


In order to assess the suitability of different modelling approaches and RANS turbulence models, it 
is worthwhile comparing the predictions to higher fidelity LES or hybrid LES results. However, the 
large thermal mass associated with the metal structure means that the mesh size and timescales 
associated with a higher fidelity simulation would be unfeasible. Therefore, the TALL-3D study has 
been split into two parts: 


Benchmark data (fluid only): A WMLES model of just the LBE fluid inside the 3D test section 
has been developed and solved (Section 5). This provides a set of benchmark data with 
well-defined boundary conditions. This WMLES data has been used to investigate and as- 
sess different modelling approaches (2D axisymmetric, 3D steady and 3D unsteady), RANS 
turbulence models (realizable k-¢ and k-w SST) and turbulent Prandtl number modifica- 
tions, as well as the sensitivity to mesh density. 


Validation exercise (full CHT): Once the most appropriate modelling approach and turbulence 
model has been confirmed for each scenario (forced and mixed convection), the model was 
extended to include the stainless steel vessel and thermal insulation (Section 6). This full 
CHT model enabled the CFD predictions to be validated against the measured data to quan- 


18 of 61 


4.2 


Planning the Analysis 


tify the accuracy of the predictions and sensitivity of the results to the material properties, 
initial and inlet boundary conditions. 


The TG03.S301.04 TALL-3D test was run at KTH in January 2017 as part of the SESAME project 
(Grishchenko et al., 2017). It is a forced to natural circulation transient with a constant electric 
power in the main heater and 3D test section. This test has been selected for the validation exercise 
because it is an open benchmark test for validation of coupled system and CFD codes with well- 
defined statistically steady flow conditions at the beginning and end of the transient. The data for 
this test was gratefully received from KTH in order to undertake this case study. 


The results of importance for the TALL-3D test case are the solid and fluid temperatures within the 
3D test section, which can be compared to the experimental measurements. This should consider 
both the temperature profile and absolute values in order to understand and assess the level of 
validation that has been achieved and the confidence that it provides for future analyses. The level 
of agreement required will depend on how this validation exercise will be used (i.e. whether it would 
be used in a chain of supporting evidence for a safety argument of a reactor system), which will 
also impact the level of quality assurance required. 


Quality Assurance 


A planned approach to quality assurance is recommended by Volume 1 (Section 4.5) in order to 
identify problems early and minimise the amount of re-work. 


An overview of the quality plan for this case study is given below: 


A suitably qualified and experienced person is selected to lead the verification of the case 
study, independently from the originator. 


The geometry of the model was verified against the details in Grishchenko et a/. (2017). 


The generation of the mesh was checked for calculated sizing, quality and a general review 
of the expected flow features. 


The boundary conditions and material properties were checked against the supplied test 
data and that they were applied correctly in the model. 


The settings for the physical and numerical models were checked to make sure they are cor- 
rectly applied in accordance with the analysis software recommendations and best practice. 


The convergence, solution behaviour and results were checked to confirm the expected flow 
phenomena are adequately modelled and resolved. 


The post-processing scripts and reported figures/numbers were checked to ensure that the 
data has been plotted correctly. 
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4.3 


Planning the Analysis oe, 


Material Properties 


The LBE properties have been taken from NSC (2015), a reference that has compiled the main 
thermophysical properties of molten LBE reported in open literature and proposed recommenda- 
tions based on the ‘best fit’ to the available data. Further information on liquid metal properties and 
additional sources of data are provided in Volume 5 (Section 2.1). The following correlations for the 
LBE fluid properties have been implemented in a UDF, where T is the temperature in Kelvin. 


Density (kg/m°): 
p = 11065.0 — 1.293T 


Dynamic viscosity (kg/m s): 
4.1 
pe = 4.94 x 10-* exp (=) 


Thermal conductivity (W/m kK): 
k = 3.284 + 1.617 x 10-77 — 2.305 x 10-°T? 
Specific heat capacity (J/kg kK): 


Cp = 164.8 — 3.94 x 10-?T + 1.25 x 10-°T? — 4.56 x 10°T~? 


The properties of American Iron and Steel Institute (AISI) 316L stainless steel have been taken 
from Incropera et al. (2011) with constant density and temperature varying properties for thermal 
conductivity and specific heat capacity. Constant density and specific heat capacity and tempera- 
ture varying thermal conductivity properties for the two types of thermal insulation (ISOVER Tape- 
Lock 7300 and Nano T Ultra) have been taken from Grishchenko et al. (2017). 
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Benchmark Data (Fluid Only) 


The fluid only CFD model is a simplified representation of the actual 3D test section. This enables 
the flow phenomena within the 3D test section to be simulated without the complexity of conju- 
gate heat transfer and the uncertainty associated with the thermal insulation and external heat 
loss. Therefore, this provides a more controlled benchmark to test, compare and assess different 
solution methodologies. 


Figure 5.1 shows the geometry of the fluid only model, which includes an adiabatic inlet and outlet 
pipe section with the 3D test section cylinder and circular inner plate. The conjugate heat transfer 
through the inner plate has been included in the model as this will affect the temperature variation 
within the 3D test section. 


—— Heater flux 


— External heat flux 


27.8 mm 
diameter 
outlet 


200 mm diameter 


200 mm 


300 mm diameter 


17mm 
diameter 
inlet 


inlet 
pipe 


Figure 5.1: Fluid only model geometry. 
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Benchmark Data (Fluid Only) 


Table 5.1 provides a summary of the inlet and boundary conditions for the mixed and forced con- 
vection cases. This includes the heat flux through the cylindrical heater and external heat loss from 
the other surfaces of the 3D test section. The external heat loss has been estimated based on the 
fluid temperature drop from the inlet to outlet in a TALL-3D test that was run with the 3D test section 
heater turned off (TG0O2.06.01). However, it is worth noting that the mixed convection external heat 
loss is much lower than the calculation based on the validation test data (Section 6). 


TG03.S301.04 Mixed convection Forced convection 
Mass flow rate (kg/s) 0.266 1.30 

Inlet temperature (°C) 201.0 235.0 
Heater power (W) 4,000 4,000 
External heat loss (W) 62.8 307.0 
Heater heat flux (W/m?) 36,396.48 35,654.31 
External heat flux (W/m?) -190.86 -933.04 


Table 5.1: Nominal flow and boundary conditions for fluid only simulations. 


Table 5.2 summarises the LBE material properties at the inlet using the correlations detailed in 
Section 4.3, which includes calculated values of the relevant non-dimensional properties discussed 
in Volume 3. The Rayleigh number (Ra = gL3(p1 — p2)Pr/vzp1) was estimated based on the 
internal length of the 3D test section to the CIP (L = 0.17 m) using the inlet temperature for p; and 
the measured ILW temperature at the CIP location for p2 (290°C and 260°C for the mixed and 
forced convection cases respectively). 


As discussed in Volume 3 (Section 2.2.1), Gr/Re” represents the ratio of buoyancy forces to mo- 
mentum forces, when n = 2 the ratio is often referred to as the Richardson number. Since this 
calculation is just intended to provide an indication of the relative strength of the buoyancy force, 
other values of n were not considered in this case. 


TG03.S301.04 Mixed convection Forced convection 
Density (kg/m?) 10451 10408 
Viscosity (kg/m s) 0.00242 0.00218 
Thermal conductivity (W/m K) 10.43 10.91 
Specific heat capacity (J/kg K) 146.9 146.2 
Prandtl number 0.0341 0.0292 

Inlet Reynolds number 8,221 44,687 

Inlet Péclet number 281 1,306 
Approx. Rayleigh number (x 10°) 337 100 
Approx. Richardson number 146 1.7 


Table 5.2: Inlet properties and non-dimensional parameters for fluid only simulations. 


The inlet Reynolds number is high in the forced convection case (Section 5.5), while the Richard- 
son number is low, which confirms that the heat transfer is expected to be dominated by forced 
convection. The Richardson number is high in the mixed convection case (Section 5.4), which sug- 
gests that the heat transfer in the 3D test section is dominated by buoyancy and natural convection, 
while the inlet Reynolds number is turbulent and so the heat transfer in the inlet pipe is driven by 
forced convection. 
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5.1 


Benchmark Data (Fluid Only) 


The fluid only cases have been solved in ANSYS Fluent using standard RANS turbulence models 
with and without modifications to the turbulent energy Prandtl number. 2D axisymmetric and 3D 
RANS simulations have been compared and assessed against an equivalent WMLES model to 
determine the most appropriate modelling approach. 


Before undertaking the RANS simulations, the predicted inlet pipe flow profiles were compared 
(Section 5.1) and the meshing approach was assessed (Section 5.2). This includes a mesh sensi- 
tivity study to confirm that the mesh is appropriate and has minimal impact on the results. 


Inlet Profile 


The inlet profile to the 3D test section was generated by extruding the mesh on the inlet face a 
small distance upstream. This separate model was solved using a periodic boundary condition 
with a specified mass flow rate at a constant temperature. The converged inlet profile was then 
imported into the 3D test section model. 


Figure 5.2 shows the velocity and turbulence profile across the inlet pipe for the RANS turbulence 
models being investigated. This has been plotted using the non-dimensional distance across the 
pipe (r/R), where r is the radial coordinate and R is the internal pipe radius. 
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Figure 5.2: Inlet profiles for fluid only simulations. 


The peak velocity predictions are generally consistent between turbulence models, with the ex- 
ception of the realizable k -¢ model of the mixed convection case, which predicts a peak velocity 
6% lower than the other two turbulence models. The prediction of turbulent viscosity ratio varies 
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significantly between turbulence models for both mixed and forced convection. The shape of the 
turbulent velocity profile is consistent with the 1/7 power law velocity profile for turbulent flow in a 
pipe (Pope, 2000). 


For the WMLES model, the steady-state RANS inlet profile was used with a random vortex method 
to generate a time varying inlet velocity fluctuation. This adds a perturbation to the mean velocity 
profile via a fluctuating vorticity field (1,000 vortices were specified). 


Meshing Approach 


The meshing approach was developed and tested on the 2D axisymmetric model before sweeping 
the 2D plane into a full 360°model with the circumferential spacing set to give approximately unit 
aspect ratio cells. A structured quadrilateral mesh was used in the 2D axisymmetric model with a 
boundary layer mesh (Figure 5.3). 


In order to maintain a consistent aspect ratio and minimise the cell count, an unstructured quadri- 
lateral mesh was generated in the central core of the inlet/outlet pipe and mid-part of the 3D test 
section (0.05m <r < 0.095m), which was then swept axially through the 3D test section. This 
means that the 3D model consists of a combination of structured and unstructured hexahedral 
cells. However, if a fully structured mesh was used, then it is likely there would be either too many 
cells at the centre or too few cells at the outer radius. 


¢ The baseline mixed convection mesh (RANS and WMLES) has a spacing of 0.3mm in the 
core jet region, which increases to 1.0mm near the ILW. The FCH of the boundary layer 
mesh was adjusted to achieve y* < 1 on all walls. 


The baseline forced convection RANS mesh is slightly more refined than the mixed convec- 
tion mesh and uses a square boundary layer mesh at the inlet pipe junction, although this 
required the FCH to increase slightly (yv* < 5). This improved the mesh resolution in the 
shear layer as the jet enters the 3D test section by reducing the radial cell spacing and elim- 
inating the skewed cells at the junction. This was more important in the forced convection 
case as the velocity gradient across the shear layer, and hence mixing, is higher. 


The forced convection WMLES mesh was created by adapting the cells in a 25mm cylinder 
up to the CIP (i.e. each cell is split into eight) to further resolve the jet shear layer. 


The coarser meshes have been created by doubling the edge spacing each time, but maintaining 
the same FCH and adjusting the boundary layer growth accordingly. This was done manually within 
the mesh generator by halving the number of cells on each edge. The finer meshes have been 
generated by adapting the baseline (and fine) mesh, which splits each cell into four (i.e. halves the 
FCH). The metrics for the baseline and coarser/finer meshes are detailed in Table 5.3, where: 


Mesh cell count 
Baseline cell count 


Mesh spacing factor = / 


Figure 5.4 plots the normalised average temperature change (%) compared to the finest mesh on 
the instrumented surfaces against mesh spacing factor for the 2D axisymmetric k-w SST turbu- 
lence model. This demonstrates that the baseline results have a low dependence on mesh density, 
and confirms that the resolution of the baseline meshes is appropriate. 
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b vane 


PAAR 


a: Mixed convection mesh: side view ction mesh: plan view 


c: Forced convection mesh: side view of inlet pipe 


Figure 5.3: Baseline mesh for fluid only simulations. 


Case Coarsest Coarse Baseline Fine Finest 
Mixed 2D mesh (cells) 7 x 103 22x10? 76x10? 303x10%? 1.2 x 10° 
Mixed mesh spacing factor 0.31 0.54 1.0 2.0 4.0 
Mixed 3D/WMLES mesh (cells) 23 x 10° 

Forced 2D mesh (cells) 11 x 10° 28 x 10° 86x 10? 345x 10% 1.4 x 10° 
Forced mesh spacing factor 0.36 0.57 1.0 2.0 4.0 
Forced 3D mesh (cells) 35 x 10° 

Forced WMLES mesh (cells) 110 x 10° 


Table 5.3: Mesh metrics for Fluid only simulations. 
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Figure 5.4: Mesh sensitivity results. 


LES mesh checks: The LES mesh resolution was reviewed based on the 3D steady RANS so- 


lution by calculating the integral length scale and Taylor microscale values (Volume 2, Sec- 
tion 3.4.2.4). Typically, 10 to 20 mesh cells are needed across /g to resolve 80% of the 
turbulent kinetic energy, so the suggested cell size is typically taken to be A = max(A, /o/10) 
for full LES resolution (further detail is available in Addad et al., 2008). 

k3/2 ki/2 


Integral length scale (m): /) = —— = ——, where C, = 0.09 (Menter, 2015) 
€ Cw 


|. k : : oe : 
Taylor microscale (m): > ~ 4/10v—, where v is the kinematic viscosity 
& 


Cell length (m): A = V/3 for modest aspect ratios, where V..) is the cell volume 


cell 


This can be checked by plotting /o/A < 10 and \/A < 1 to identify areas that are under- 
resolved. The integral length scale and Taylor microscale values were compared for the 
mixed and forced convection cases and the worst-case values are plotted in Figure 5.5 with 
increased limits for clarity. This shows that the mesh size should be sufficient to resolve over 
80% of the turbulent kinetic energy in the bulk of the test section for both forced and mixed 
convection cases. The LES resolution quality should be checked using the LES results to 
assess the proportion of turbulent kinetic energy that was resolved. 


2 
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a: Mixed convection, integral length scale / A b: Forced convection, Taylor microscale / A 


Figure 5.5: RANS estimate of LES resolution on symmetry plane. 
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Thermal Mass of Inner Plate 


Although it is more accurate to include the conjugate heat transfer through the inner plate, it is im- 
portant to understand the impact of the thermal mass (mc,) in the solid on the timescales required 
to achieve a steady-state solution’. This is discussed in Volume 2 (Section 2.1.1.2). 


In order to understand the temperature gradients through the inner plate and determine the con- 
duction timescales, it is necessary to calculate the Biot number (Bi = hAL/k,) and Fourier number 
(Fo = a,t/L?) for the inner plate respectively, where his the fluid Heat Transfer Coefficient (HTC), 
k, is the solid thermal conductivity, L is a representative length (plate volume / area), t is the 
elapsed time and a, is the solid thermal diffusivity (k/pcp). These values can be easily extracted 
from the CFD model, except for the HTC, which needs to be calculated for the upper and lower 
surface of the plate using one of the methods described in Volume 2 (Section 3.4.6). 


Since there are large temperature variations in the fluid below the inner plate due to the jet, there is 
no representative fluid temperature. Therefore, the thermal superposition method has been used 
to predict the HTC on the plate surface. This was achieved by solving a 2D axisymmetric model 
without the plate in order to calculate the adiabatic surface temperature, and then re-solving the 
CFD model with a 1 °C temperature change to the upper and lower surface separately. This enabled 
the HTC, Biot number (B/) and Fourier number (Fo) to be calculated (Table 5.4). 


Mixed convection Forced convection 
Average HTC (W/m? kK) 800 4000 
Biot number 0.2 1.06 
Fourier number (t= 1s) 0.193 0.189 
Time period, to.s (Ss) 17.6 3.5 


Table 5.4: Unsteady conduction through circular inner plate. 


Since the Biot number is larger than 0.1, there will be a thermal gradient through the inner plate 
and a ‘lumped capacitance’ approach is not valid. However, an indication of the time period of the 
plate can be made using this method, as follows: 


0 _ T~Teo _ Bio _ .-(wE)t pV Cp 


pVcp rs to.5 _——— 


ae ri : hs 


In(0.5) 


If the plate starts at 10°C above the surrounding fluid, to.5 is the time taken for the plate to reduce 
to 5°C above the surrounding fluid (i.e. the difference is halved). This suggests that the time taken 
for the plate to reach steady-state conditions is significant when compared to the time step required 
for unsteady CFD models. However, since the unsteady simulations start from a converged steady 
solution, the impact on the temperature predictions is expected to be acceptably small. 


Due to the large thermal mass and long timescales associated with the inner plate, care needs to 
be taken with the unsteady simulations to ensure that the inner plate temperatures have properly 
converged. There is a risk that the plate temperature could slowly heat up or cool down as the 
solution progresses in time. Therefore, the plate needs to be monitored during the solution to 
confirm that the temperatures have stabilised and are fluctuating about a constant mean value. 


It is worth noting that for liquid metal the thermal mass of the fluid in the 3D test section is significantly higher than the inner 
plate. However, since the fluid is moving the residence time is a more relevant factor for solution timescales. 
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Mixed Convection 


The 2D axisymmetric and 3D steady k-w SST models have been solved using the baseline mesh 
with gravity applied vertically downwards. These models were initialised with zero velocity and low 
turbulence, and solved using the SIMPLE pressure-velocity coupling scheme with second-order 
discretisation and double precision. The under-relaxation parameters were set to the default values 
except for the energy equation, which was reduced to 0.99 in order to stabilise the solution. 


The temperature was monitored at a number of thermocouple locations during the solution, and 
convergence was considered acceptable when the monitors and scaled residuals had levelled off. 


Figure 5.6 shows velocity and temperature contours within the 3D test section, which demonstrate 
that the flow in the 3D steady model is circumferentially uniform. These are plotted on a symmetry 
plane through the centre of the model, and on a horizontal plane four inlet pipe diameters (4D) 
above the base plate, where D = 0.017 m. The inner plate is located 10D above the base plate. 


A comparison of the temperatures on the instrumented surfaces confirms that the 2D axisymmet- 
ric and 3D steady results are within 0.01 K of each other, and so a 2D axisymmetric model is 
appropriate for this case. 
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Figure 5.6: Mixed convection: 3D steady k-w SST results. 


The WMLES model has been solved using a time step of 0.0012 s to ensure that the CFL number 
is less than 1. The temperature was monitored at a number of locations within the 3D test section 
to confirm that the RANS solution had washed through and the fluctuations had stabilised. The 
time-averaged mean velocity and temperature data were generated by solving for a further 20,000 
time steps, and the mean values were then circumferentially averaged to create a single set of data 
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for each surface extracted. 


The instantaneous velocity and temperature contours on the symmetry plane (Figure 5.7) show the 
unsteadiness in the low velocity jet with little fluctuation in the thermally stratified region, while the 
time-averaged results show that the mean jet profile is symmetric and consistent with the steady 
RANS results (Figure 5.6). 
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Figure 5.7: Mixed convection: WMLES results. 


Jet development: Figure 5.8 compares the RANS velocity and temperature profiles with and 
without the Kays modification for the jet on the 2D, 4D, 6D and 8D planes above the base plate 
against the equivalent circumferential and time-averaged WMLES results. The radial coordinate 
has been normalised by the inlet pipe diameter (r/D). 


The small notches in the WMLES results at r/D = 2 are due to the change in the mesh from 
structured to unstructured hexahedrals, and are due to the number of cells selected at each radial 
coordinate as part of the circumferential averaging routine. 


The results show that the velocity profile is generally captured by the RANS turbulence models, 
except for the realizable k - ¢ model, which decays too quickly. In general, the temperature profile 
is predicted by all models, and the Kays modification to the turbulent energy Prandtl number has 
a relatively small impact on the results. The WMLES model predicts the negatively buoyant jet 
reaches slightly higher than the RANS predictions, as shown in the 8D plane temperature profiles. 
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Figure 5.8: Mixed convection: Jet velocity and temperature profiles. 
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Surface temperatures: Figure 5.9 compares the RANS temperature profiles with and without 
the Kays modification on the instrumented surfaces (BP, CIP, IPT and ILW) of the 3D test section 
against the equivalent circumferential and time-averaged WMLES results. 


The results show that the 2D axisymmetric RANS turbulence models predict the correct trends in 
the surface temperatures, although there is some variation between the different models. 


Overall, the k - w SST with full buoyancy terms and without the Kays modification provides the clos- 
est match to the WMLES results (within 1 to 2°C), as buoyancy effects are significant in the mixed 
convection case. The realizable k -€ model includes buoyancy effects on the turbulence produc- 
tion (k), but not the turbulence dissipation (¢) by default (Volume 3, Section 2.2.4.3). Although the 
realizable k-«€ wall temperature predictions are similar to the k-w SST with full buoyancy terms, 
the jet development, discussed below, is not as accurately predicted. 


Although the Reynolds Stress Model (RSM) turbulence model solves the transport equations for 
each of the Reynolds stresses directly, the comparison to the WMLES results is not particularly 
good especially on the IPT and ILW surfaces. Therefore, the results suggest that it is not appropri- 
ate in this case. 
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Figure 5.9: Mixed convection: Instrumented surface temperature. 
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Forced Convection 


The 2D axisymmetric model has been solved steady-state using the realizable k-¢ and k-w SST 
turbulence models with gravity applied vertically downwards. These models were initialised with 
zero velocity and low turbulence, and solved using the SIMPLE pressure-velocity coupling scheme 
with second-order discretisation and double precision. The under-relaxation parameters were set 
to the default values including the energy equation. 


A comparison of the velocity and temperatures within the 3D test section (Figure 5.10) show that 
there is a large difference in the flow profile, and hence temperatures, between these two RANS 
turbulence models due to the change in the flow recirculation and streamlines. 
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Figure 5.10: Forced convection: 2D axisymmetric results. 
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The 3D model has been solved using both a steady RANS and Unsteady Reynolds-Averaged 
Navier-Stokes (URANS) approach using the realizable k - ¢ and k-w SST turbulence models. The 
velocity and temperature contours (Figure 5.11) show that the jet is asymmetric in the 3D steady 
solution, which reduces the amount of recirculation in the 3D test section and increases the bulk 
fluid temperature compared to the 2D axisymmetric results. The URANS solution (time step = 
0.004 s) simulates the unsteady fluctuating motion of the jet, and confirms that the jet is not steady. 
However, when the URANS predictions are averaged over time, the time-averaged profile of the jet 
is symmetric with a similar flow pattern to the equivalent 2D axisymmetric streamlines (Figure 5.10). 


b: Steady RANS temperature (°C) 


c: Time-averaged URANS velocity magnitude (m/s) d: Time-averaged URANS temperature (°C) 


Figure 5.11: Forced convection: 3D realizable k - € results. 


The WMLES model has been solved using a time step of 0.0002 s to ensure that the CFL number 
is less than 1. The temperature was monitored at a number of locations within the test section to 
confirm that the RANS solution had washed through and the fluctuations had stabilised. The mean 
velocity and temperature data were generated by solving for a further 20,000 time steps, and the 
mean values were then circumferentially averaged to create a single set of data for each surface 
extracted. 


The instantaneous velocity and temperature contours on the symmetry plane (Figure 5.12) show 
the unsteadiness in the high velocity jet and recirculating region within the 3D test section, while 
the time-averaged results show that the mean jet profile is symmetric and consistent with the time- 
averaged URANS results (Figure 5.11). 
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Figure 5.12: Forced convection: WMLES results. 


Jet development: Figure 5.13 compares the 2D axisymmetric, 3D RANS and 3D URANS tem- 
perature profiles with and without the Kays modification for the jet on the 2D, 4D, 6D and 8D planes 
above the base plate against the equivalent circumferential and time-averaged WMLES results. 
The radial coordinate has been normalised by the inlet pipe diameter (r/D). 


The results show that the velocity and temperature profile is generally captured by the RANS tur- 
bulence models, and the Kays modification to the turbulent energy Prandtl number has a relatively 
small impact on the results. The 2D axisymmetric results predict less mixing in the jet, which leads 
to increased velocities and reduced temperatures in the core. 


Surface temperatures: Figure 5.14 compares the 2D axisymmetric, 3D RANS and 3D URANS 
temperature profiles with and without the Kays modification on the instrumented surfaces (BP, CIP, 
IPT and ILW) of the 3D test section against the circumferential and time-averaged WMLES results. 


The results show that the ILW temperature profile for the k - w SST turbulence model is significantly 
different to the WMLES due to the change in the recirculation structure and location of the separa- 
tion point (Figure 5.10). This suggests that the realizable k - ¢ model is a more appropriate RANS 
turbulence model for the forced convection case. 


Although the 2D axisymmetric results are consistent with the 3D RANS and WMLES predictions, 
the temperatures are approximately 2°C lower due to the increased recirculation below the CIP. 
This suggests that a 2D axisymmetric approach is not appropriate in this case and could lead to 
an under-prediction of key temperatures. 


The 3D RANS and URANS results are generally similar to each other, which suggests that, in this 
case, the steady RANS solution is approximately equivalent to a snapshot in the unsteady flow. 
However, since the steady solution temperatures vary with each iteration it is hard to accurately 
compare different predictions, and such variation does not have physical significance (it is a feature 
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Figure 5.13: Forced convection: Jet velocity and temperature profiles. 
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of the numerical solution algorithm), so a URANS solution approach is more appropriate for this 
case due to the flow instability. 


The Kays modification has a relatively small impact on the temperature predictions, but generally 
improves the agreement with the WMLES results. In some areas the results without the Kays 
modification are closer to the WMLES results, such as the ILW surface, and so it is worthwhile 
assessing the sensitivity of the results to the turbulent Prandtl number. In the future, it may be 
worthwhile investigating whether an AHFM approach to model the turbulent heat flux improves the 
temperature predictions. 


Overall, the realizable « -¢ URANS model with the Kays modification provides the closest match to 
the WMLES results (< 1 °C) for the forced convection case. 
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Figure 5.14: Forced convection: Instrumented surface temperature. 
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Using the modelling approach identified as providing the best results for each scenario in the 
benchmark comparisons (Section 5), the model was extended to include the stainless steel vessel 
and thermal insulation. This full CHT model enables the CFD predictions to be validated against the 
TALL-3D thermocouple measurements to quantify the accuracy of the predictions and sensitivity 
of the results to the material properties, inlet and boundary conditions. 


Figure 6.1 shows the geometry of the full CHT model (Grishchenko et a/., 2017), which includes 
the stainless steel inlet/outlet pipe, test section components and heater, as well as the surrounding 
insulation. In addition, the model has been extended upstream/downstream to the thermocouple 
measurement locations and includes the air cavities between the insulation and stainless steel. 


The only difference in the 3D test section geometry compared to the fluid only model is that the 
heater length is 120mm (instead of 116mm); KTH confirmed that one of the 3D test section rope 
heaters was replaced in 2016 before the SESAME test campaign. 
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Figure 6.1: Full CHT model geometry. 


The temperature varying fluid (LBE), stainless steel and insulation material properties are de- 
tailed in Section 4.3, while the temperature varying air properties were taken from Incropera et al. 
(2011). Table 6.1 provides a summary of the inlet and boundary conditions for the mixed and forced 
convection cases. This includes the heat flux through the cylindrical heater, external ambient air 
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temperature, outlet temperature and calculated external heat loss (Qjoss)- 


Tr 
ue Cp,out + Cp,in 
Qioss = Q heater a w | cpdT ~~ Q heater W ( P 7 Pr ) (Tait —_ T; ) 
Tin 


This approximation is appropriate because the specific heat capacity is effectively linear over the 
temperature range within the 3D test section (Rogers and Mayhew, 1992)'. 


The outlet temperature includes the offset calculated in Table 3-13 of Grishchenko ef al. (2017), 
which provides a correction to the differential reading between the inlet and outlet thermocouples. 
This was calculated by measuring the temperature difference across the 3D test section during 
forward and reverse flow conditions by assuming that the thermal losses are the same in both 
directions. This should reduce the uncertainty of the differential temperature measurement across 
the 3D test section. 


TG03.S301.04 Mixed convection Forced convection 
Mass flow rate (kg/s) 0.266 1.31 

Inlet temperature (°C) 201.82 235.41 
Heater power (W) 4,025.86 4,048.44 
External bottom air temperature (°C) 23.17 22.82 
External middle air temperature (°C) 28.12 28.36 
Outlet temperature (°C) 295.52 254.83 
Calculated external heat loss (W) 390.7 333.7 


Table 6.1: Nominal flow and boundary conditions for full CHT simulations. 


Meshing Approach 


The mesh in the 3D test section was closely matched to the fluid only simulation mesh (Sec- 
tion 5.2), with the additional surrounding solid components meshed using structured/unstructured 
quadrilateral/hexahedral cells for the 2D axisymmetric mixed convection and 3D forced convection 
cases respectively (Figure 6.2). The fluid and solid regions in the forced convection case are gen- 
erally conformal with a small number of non-conformal interfaces used in specific regions, such as 
the centre of the CIP. The total number of cells and LBE fluid cells are summarised in Table 6.2 for 
each case. 


It is worth noting that the solid mesh is probably more refined than is necessary in order to resolve 
the thermal gradients in the solid. This eliminates the need to undertake a mesh sensitivity study 
and is acceptable for this validation study, especially for the short solution times of the 2D axisym- 
metric models. However, significant mesh savings could be achieved by using a non-conformal 
fluid-solid mesh interface and a coarser solid mesh. This would be more appropriate for practi- 
cal engineering applications, but the mesh resolution should be checked using a mesh sensitivity 
study. 


It is worth noting that simpler expressions can result in significant errors in the calculated heat loss. For the natural convec- 
tion case, Qheater — Wp, in( Tout — Tin) is 6 % lower than the actual heat loss, while Qheater — W(Cp,out Tout — Cp,in Tin) is 
66 % higher than the actual heat loss. 
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b: Top of 3D test section 


c: Bottom of 3D test section 
a: Outer wall of 3D test section 


Figure 6.2: Full CHT solid mesh. 


Case Total LBE fluid 
2D mixed convection mesh (cells) 174 x 10° 85 x 10° 
3D forced convection mesh (cells) 68 x 10° 38 x 10° 


Table 6.2: Mesh metrics for full CHT simulations. 


Solution Approach 


The full CHT cases have been solved in ANSYS Fluent using the following RANS turbulence mod- 
els based on the results of the fluid only simulations (Section 5): 


* Mixed convection (2D axisymmetric, steady): k-w SST turbulence model with full buoyancy 
terms and default turbulent energy Prandtl number. 


« Forced convection (3D, unsteady): Realizable & - ¢ turbulence model with Kays modification 
to the turbulent energy Prandtl number. 


The inlet profile was generated using the same periodic section approach as the fluid only model 
(Section 5.1) by extruding the mesh on the inlet face a small distance upstream. This separate 
model was solved with a specified mass flow rate and a constant heat flux applied to the inlet pipe 
wall of -517.24 W/m? (mixed) and -659.18 W/m? (forced) convection cases. This was calculated 
from the full CHT CFD model based on the heat loss through the pipe wall and insulation to the am- 
bient air with the HTC and ambient temperature estimated using the correlations in Section 6.2.1. 
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The converged inlet profile was then imported into the 3D test section model. 


The temperatures in the solid and fluid were solved together in a single solver using a coupled 
monolithic approach, as detailed in Volume 2 (Convection, Radiation and Conjugate Heat Trans- 
fer) with full contact assumed between different solid components. Gravity was applied vertically 
downwards, and the solution approach was consistent with the fluid only simulations (Section 5). 


The air cavity regions were solved using the Boussinesq approximation with density based on the 
average temperature in each cavity. This was necessary to enable the model to solve multiple 
separate sealed fluid regions and ensure that the pressure/velocity variations remained stable. Ra- 
diative heat transfer across the air cavities was included in the model using the Discrete Ordinates 
(DO) radiation model with surface emissivity of 1 to maximise radiative heat transfer, and hence 
external heat loss. 


The thermal conductivity of the NANO T Ultra insulation on the top of the 3D test section was 
increased to 0.155 W/mK to take into account the thermocouple tubes that pass through the 
insulation (Figure 3.2). This has been estimated based on the proportion of stainless steel present. 


External Heat Transfer 


The heat loss from the external surface of the insulation is simulated in the CFD model by apply- 
ing a convection (HTC and ambient temperature) and thermal radiation (emissivity and radiation 
temperature) boundary condition: 


¢ The insulation is coated in a shiny reflective coating, so the external emissivity was estimated 
as 0.33. 


* The external ambient temperature was set to the measured external bottom air temperature 
for the bottom and side surfaces, and the external middle air temperature for the top surfaces. 


« The HTCs were initially calculated using empirical correlations for natural convection on 
heated horizontal and vertical plates (Incropera et al., 2011): 


Lower surface of heated horizontal plate, Nu = 0.52 Ra’/® 


Upper surface of heated horizontal plate, Nu = 0.15 Ra?/? 


2 
0.387 Rai/® 


4/9 
(+ 88°" 


Vertical heated flat plate, Nu = | 0.825 4 


ilm : w !oo Tw Too N Kim 
Ra = 2 Prim LT — Teo) where Thm = ret Feo) ~ ) an nTce MOM 


OLFilm Vtilm L 
The external boundaries were split into the following surfaces: Horizontal plate (bottom surface, 
outer and inner top surface split at outer radius of heater) and vertical plate (inlet pipe, outer 
surface and outlet pipe). 


Some iteration was required in order to ensure that the predicted wall temperatures and calcu- 
lated HTCs are consistent. However, since the external surfaces are complex and do not match 
the idealised conditions associated with the correlations, these values can only be considered as 
estimated values. 
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Mixed Convection 


This section describes the results of the full CHT mixed convection results. A baseline case was 
initially developed and solved using the empirically calculated boundary conditions detailed in Sec- 
tion 6.2. A sensitivity analysis was then conducted to understand the impact of material properties, 
inlet and boundary conditions on the temperature predictions (Section 6.3.2) before generating a 
best estimate solution of the mixed convection case (Section 6.3.3). 


Baseline Results 


The external boundary conditions calculated using the empirical correlations detailed in Sec- 
tion 6.2.1 and applied to the baseline case are listed in Table 6.3. 


Baseline case Inlet pipe Base Side Topouter Topinner Outlet pipe 
External HTC (W/m? kK) 3.0 1.2 a1 3.6 6.4 3.6 
Air Temperature (°C) 23.17 23.17 23.17 28.12 28.12 28.12 


Table 6.3: External boundary conditions for baseline model. 


The velocity and temperature contours are plotted in Figure 6.3, which demonstrates that the re- 
sults are consistent with the fluid only simulation (Figure 5.6). As expected, the temperature drops 
to just above ambient through the steel casing and insulation with the largest thermal gradient 
predicted in the insulation due to the low thermal conductivity. 


a: Velocity magnitude (m/s) b: Temperature (°C) 


Figure 6.3: Full CHT, mixed convection results. 


The baseline results are plotted against the fluid only results and TALL-3D test data in Figure 6.4. 
The test data include the IPT thermocouple offsets in Table 3-14 of Grishchenko ef al. (2017), 
which were calibrated using a forced circulation test with no heating in the 3D test section (i.e. fully 
mixed constant temperature) to reduce measurement uncertainty. 


The outlet temperature is 9.5°C higher than the measured value as the predicted external heat 
loss (146 W) is lower than expected. This is highlighted by the higher than measured CIP and 
upper IPT/ILW temperature predictions. 


The difference between these predictions and the test data are larger than expected, especially 
given the close agreement between the equivalent RANS and WMLES fluid only results. It therefore 
suggests that some aspects of the baseline case do not appropriately represent the actual test 
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conditions, and so a sensitivity analysis has been undertaken (Section 6.3.2) to understand which 
properties make the most difference and whether a credible variation in them can improve the 
agreement. 
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Figure 6.4: Full CHT, mixed convection: Comparison of baseline results and test data. 


Sensitivity Analysis 


As with any experimental test, there are a large number of uncertainties in the experimental mea- 
surements and test conditions that can significantly impact the accuracy of simulation results. The 
known uncertainties in the experimental measurements are detailed in Grishchenko et a/. (2017) 
and a detailed Verification, Validation and Uncertainty Quantification (VVUQ) exercise has been 
undertaken by KTH (Jeltsov et a/., 2019). 


Since this case study is specifically focused on validation against known test data, rather than a 
blind benchmark, a formal Uncertainty Quantification (UQ) exercise has not been conducted and a 
Sensitivity Analysis (SA) has been used to investigate and identify the key parameters that impact 
the temperature predictions. This was achieved by varying parameters one at a time using values 
at the upper and lower end of the plausible range, and then combining the more influential inputs 
to provide a best estimate solution. 


This kind of parameter study is suitable where the physical significance and interpretation of the 
varied inputs and the resulting output is straightforward to understand; other SA and UQ techniques 
are discussed in more detail in Volume 4 (Confidence and Uncertainty), and an example of their 
automation is described in Study B (Fuel Assembly CFD and UQ for a Molten Salt Reactor). 
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Since both SA and UQ methods require a large number of cases to be run, the mixed convection 
case (2D axisymmetric model) is more suitable for undertaking a sensitivity analysis than the forced 
convection case (3D URANS model) due to the significantly shorter solution timescales associated 
with it. 


The main sources of uncertainty in the model that were considered during the sensitivity analysis 
are discussed below with the upper and lower bounds that were tested in the model together with 
their justification listed in Table 6.4. It is worth noting that this is not considered a comprehensive 
list of all of the uncertainties in the model, and does not simultaneously cover variations related to 
the modelling itself. 


LBE material properties: As described in Section 4.3, NSC (2015) proposed correlations for the 
main thermophysical properties of molten LBE based on the ‘best fit’ to the available data. 
In addition, the maximum percentage discrepancy from the ‘best fit’ correlation has been 
quantified based on the available data and used as bounds for the sensitivity analysis. The 
variation in some of the properties is substantial, because there are few modern available 
data sources, and the observed values can vary significantly with even small quantities of 
impurities. 

Inlet conditions: The inlet temperature and mass flow rate at the model inlet are based on a ther- 
mocouple and Coriolis flow meter measurement respectively. The measurement uncertainty 
associated with these values has been estimated by KTH in Grishchenko et al. (2017). 


Heat input: The heat flux supplied into the 3D test section depends on the effective power of 
the electric rope heaters and the contact with and conduction through the steel wall. The 
uncertainty associated with the effective power measurement was estimated in Grishchenko 
et al. (2017), and the impact of a contact resistance between the heater and steel wall has 
been considered using a thin wall to account for conduction through an air gap. The uniformity 
of the heater power and thermal conductivity of the stainless steel has not been considered. 


External heat loss: The external heat loss from the 3D test section depends on the effective 
thermal resistance of the steel wall, air cavities, insulation and external HTC. The thermal 
conductivity of the insulation can vary significantly from the manufacturer’s specification due 
to how it has been installed, the level of compression and the number of penetrations or 
air gaps that exist. The external HTC could be much larger than calculated due to forced 
convection flow over the insulation surface (caused by a fan located at the top of the TALL- 
3D facility room that was used to cool the secondary loop) or increased effective external 
surface area due to the thermocouple tube penetrations. Therefore, only increased values 
compared to the initial empirical correlations have been considered as part of the sensitivity 
analysis, estimated based on engineering judgement. In addition, the impact of replacing the 
air cavities with the insulation material local to the cavity has been considered. 


Temperature measurements: The uncertainty of each thermocouple is + 1.5°C (Grishchenko 
et al., 2017), but since each data point is an average of four thermocouples at different az- 
imuthal locations the uncertainty will be reduced. The test data uncertainty has been calcu- 
lated based on an average of the temperature range for the BP, CIP and ILW thermocouples 


as +0.85°C. Since the IPT thermocouples were calibrated during the measurement cam- 


paign and offsets calculated, their uncertainty was set to half this value (+ 0.425 °C). 
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Case Property Lower bound Upper bound Justification 
1 LBE conductivity (%) -15.0 +15.0 NSC (2015) 
2 LBE specific heat capacity (%) -5.0 +5.0 NSC (2015) 
3 LBE viscosity (%) -8.0 +8.0 NSC (2015) 
4 Inlet temperature (°C) -1.5 +1.5 Grishchenko ef al. (2017) 
5 Mass flow rate (%) -3.0 +3.0 Grishchenko et al. (2017) 
6 Heater power (W) -70 +70 Grishchenko et a/. (2017) 
7 Heater/steel air gap (mm) - 1.0 Engineering judgement 
8 Air cavities removed - Local insulation Engineering judgement 
9 ISOVER conductivity (W/m kK) 0.096 (x2) 0.192 (x4) Engineering judgement 
10 NANO conductivity (W/m kK) 0.03868 (x2) 0.07736 (x4) Engineering judgement 
11 Top inner insulation (W/m kK) 0.775 (x5) 1.55 (x10) Engineering judgement 
12 All insulation conductivity As 9,10 &11 As 9,10 & 11 Engineering judgement 
13 HTC sides (W/m? k) = 15.5 (x5) Engineering judgement 
14 HTC base (W/m? K) = 6.0 (x5) Engineering judgement 
15 HTC top outer (W/m? kK) - 17.78 (x5) Engineering judgement 
16 HTC top inner (W/m? kK) . 95.6 (x15) Engineering judgement 
17 HTC all surfaces (W/m? kK) 2 As 13, 14, 15 & 16 Engineering judgement 
18 Emissivity external/top inner - 0.8/1.0 Engineering judgement 
19 Combined high heat loss - As 12, 17 & 18 Engineering judgement 


Table 6.4: Summary of sensitivity analysis parameters considered. 


In order to assess the model response to the sensitivity analysis parameters, it is necessary to 
define quantities of interest or SRQs that can be compared in both the simulation (predicted) and 
experiment (measured). The following values have been chosen in this study to assess the sensi- 
tivity of the model to each parameter, and identify the most appropriate values that should be used 
in the best estimate analysis: 


LBE outlet temperature: This parameter validates whether the overall heat balance in the model 
is correct (i.e. the combination of heat input, heat losses and energy in the LBE fluid). How- 
ever, since this parameter incorporates multiple different factors, it cannot be relied on solely 
for validation as it could be matched by mutual compensation of different errors. 


External insulation temperature: The four thermocouples (two axial locations) on the outer sur- 
face of the insulation have been extracted from the simulation and averaged. This provides a 
single value to assess whether the external heat loss through the sides of the 3D test section 
is appropriate. 


Wall surface temperatures: The predicted temperature at each thermocouple location has been 
extracted from the simulation and averaged for each instrumented surface (BP, CIP, IPT and 
ILW). This provides a single value for each surface that can be compared to the test data to 
quantify the overall comparison to the 3D test section thermocouple measurements. 
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The sensitivity analysis results are detailed in Table 6.5 as the difference between the CFD simu- 
lation results and measured values (AT, = Tsim — Tmeas). The results demonstrate that: 

* The temperature predictions for the baseline case are significantly higher than the measure- 
ments, as the heat loss from the 3D test section is too low, and/or the fraction of heat flux 
directed from the rope heater towards the 3D test section is too high. 

¢ Large variations in the temperature predictions can be achieved due to the sources of un- 
certainty in the experiment. Therefore, the uncertainty in the modelling approach is much 
smaller than the uncertainty in the boundary conditions. 

¢ The insulation conductivity and external HTC must both be increased in order to reduce the 
overall thermal resistance sufficiently to achieve a level of heat loss that agrees with the data. 

The results of this sensitivity analysis are consistent with the results of the more detailed VVUQ 

study undertaken by Jeltsov ef a/. (2019), which also demonstrated the difficulty in modelling the 

heat losses appropriately. 
Property Outlet External BP CIP IPT ILW Average 
LBE conductivity 9.5/9.5 16.2/16.6 -0.3/1.8 6.7/7.9 2.0/4.7 0.9/3.2 5.9/7.3 
LBE specific heat capacity 9.5/9.5 16.4/16.4 0.7/0.9 7.4/7.4 3.3/3.5 1.9/2.2 6.5/6.6 
LBE viscosity 9.5/9.5 16.4/16.4 0.8/0.8 7.4/7.4 3.3/3.4 2.0/2.1 6.5/6.6 
Inlet temperature 8.0/11.0 16.2/16.7 -0.7/2.3 5.9/8.9 1.9/4.9 0.5/3.6 5.3/7.9 
Mass flow rate 12.7/6.5 16.5/16.3 0.5/1.1 10.6/4.4 4.8/2.1 3.3/0.9 8.1/5.2 
Heater power 7.6/11.3 16.3/16.5 0.6/1.0 5.6/9.1 2.5/4.3 1.2/2.9 5.6/7.5 
Heater/steel air gap 7.8 77. 2.5 3.6 3.1 1.5 15.9 
Air cavities removed 9.5 15.9 0.8 7.4 3.4 2.0 6.5 
ISOVER conductivity 8.4/6.4 26.5/39.9 -0.1/-1.6 6.6/5.3 2.7/1.5 1.4/0.4 7.6/8.6 
NANO conductivity 8.8/7.7 25.0/39.0 0.5/0.0 6.8/5.8 3.0/2.4 1.7/1.2 7.6/9.3 
Top inner insulation 7.3/6.3 16.4/16.4 0.8/0.8 6.2/5.5 2.9/2.6 1.8/1.7 5.9/5.5 
All insulation conductivity 5.5/1.2 34.9/61.5  -0.5/-2.7 4.7/1.8 1.8/-0.3 0.8/-0.9 7.9/10.9 
HTC sides 9.4 -4.3 0.7 7.3 3.3 1.9 3.0 
HTC base 9.5 16.4 0.8 7.4 3.4 2.0 6.6 
HTC top outer 9.5 16.4 0.8 7.4 3.4 2.0 6.6 
HTC top inner 9.1 16.4 0.8 7.2 3.3 2.0 6.5 
HTC all surfaces 9.0 -4.3 0.7 7.0 3.2 1.9 2.9 
Emissivity external/top 9.3 6.2 0.7 7.2 3.3 2.0 4.8 
Combined high heat loss -7.6 14.2 -4.3 -3.6 -3.0 -2.7 -1.2 
Baseline case 9.5 16.4 0.8 7.4 3.4 2.0 6.6 
Best estimate case -1.2 4.2 -0.6 -0.4 -0.5 -0.9 0.1 


Table 6.5: Difference between CFD results and test data, AT, (lower/upper bound). 


The final row in Table 6.5 shows the difference between the test data and a best estimate case, 
where several adjustments to the input parameters were made simultaneously with the objective 
of minimising error in all locations. This was done by inspecting the sensitivity analysis results, and 
manually selecting adjustments that clearly improved agreement between the CFD predictions and 
test data. The following adjustments were selected: 


« An increase in LBE mass flow rate of 3% (upper bound of Case 5). 
« An increase in all insulation conductivity (lower bound of Case 12). 
« An increase in the external HTC on all surfaces (upper bound of Case 17). 
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The CFD results for each of these adjustments is compared both individually, and in combina- 
tion with the others, in Figure 6.5 as the temperature difference compared to the baseline results 
(ATp = Tease — Thase)- There is a large change in temperature on some of the surfaces, so temper- 
ature difference plots are used to more clearly show the comparison with the test measurements. 
This demonstrates the impact of each adjustment, and confirms that the best estimate solution is 
within the measurement uncertainty for most thermocouples. 
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Figure 6.5: Full CHT, mixed convection: Sensitivity analysis temperature difference, AT,. 


Best Estimate Results 


The predicted temperatures for the best estimate case, together with the individual adjustments, 
are compared to the TALL-3D measurements in Figure 6.6. The best estimate results show that 
when the mass flow rate and external heat losses from the 3D test section are increased, the CFD 
results closely match all internal and external thermocouples. 


The best estimate results demonstrate that the changes to the boundary conditions significantly 
improves the temperature profile on the IPT and ILW surfaces, as well as the CIP. This suggests 
that the higher heat loss is more representative of the real geometry, and a 2D axisymmetric CFD 
approach can be used to simulate mixed convection within the 3D test section. 


As shown in Figure 6.5 and Table 6.5, some differences still remain between the CFD predictions 
and measurement data. This suggests that the boundary conditions for the best estimate case 
could be further improved in order to better fit the experimental data. However, it is worth noting 
that: 
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* Differences exist between the 2D axisymmetric and WMLES results (Figure 6.4). Therefore, 
any differences associated with the modelling approach should be taken into account as part 
of any further model calibration process. 


« Parameters that give the best fit to the test data are not necessarily the ‘correct values’, as 
they could be matched by mutual compensation of different errors. 


However, due to the large number of uncertainty parameters, this would require a detailed op- 
timisation study involving a significant number of simulations. This could be achieved using an 
automated model calibration process to determine the best fit inputs. This is a form of optimisation, 
which is described in Volume 4. 


The purpose of this case study is to validate the capability of CFD to simulate mixed convection 
within liquid metal, rather than to interpret the details of the non-fluid aspects of the experimental 
setup, so the level of agreement that has been achieved with the best estimate case is considered 
appropriate. 
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Figure 6.6: Full CHT, mixed convection: Best estimate results. 
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6.4 


Validation Exercise (Full CHT) CRF 


Forced Convection 


This section describes the results of the forced convection model using the baseline and best 
estimate case boundary conditions detailed in Section 6.3, and compares the results of the 2D 
axisymmetric, 3D steady RANS and 3D URANS results. 


In addition, this tests the parameters that have been chosen for the best estimate case because if 
the best estimate parameters also provide a better prediction for the forced convection case, then 
it provides more confidence that the changes are appropriate. In terms of a numerical process, 
the mixed convection case is the ‘training’ or ‘model calibration’ dataset and the forced convection 
case can be considered a validation dataset. 


The external boundary conditions calculated using the correlations in Section 6.2.1 and applied 
to the baseline case are listed in Table 6.6. The baseline and best estimate insulation thermal 
conductivity matches the mixed convection simulation, while the inlet mass flow rate has been 
increased by 0.9% (rather than 3%) in line with the measurement uncertainty estimated by KTH 
in Grishchenko et a/. (2017) for higher flow rates. 


Inletpipe Base Side Topouter Topinner Outlet pipe 


Baseline case HTC (W/m? kK) 22 1.3 25 a7 6.5 3.4 
Best estimate case HTC (W/m? k) 3.2 6.0 15.5 17.78 95.6 3.4 
Air Temperature (°C) 22.82 22.82 22.82 28.36 28.36 28.36 


Table 6.6: External boundary conditions for forced convection models. 


The velocity and temperature contours for the 3D URANS model are plotted in Figure 6.7, which 
demonstrates that the results are consistent with the fluid only simulation (Figure 5.11). As for 
the mixed convection case, the temperature drops to just above ambient through the steel casing 
and insulation with the largest thermal gradient predicted in the insulation due to the low thermal 
conductivity. 
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Figure 6.7: Full CHT, forced convection URANS results. 


The full CHT results are plotted against the fluid only results and TALL-3D test data in Figure 6.8, 
which include the IPT temperature offsets in Table 3-14 of Grishchenko et a/. (2017). Since the 
mass flow rate for the forced convection case is approximately 5 times higher than the mixed con- 
vection case, the temperature rise is 5 times smaller. This means that smaller differences between 
simulations are more significant. 
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Validation Exercise (Full CHT) 
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Figure 6.8: Full CHT, forced convection: Comparison of CFD results and test data. 


The results demonstrate the impact of different modelling approaches on the monitored surface 
temperatures: 


« The 2D axisymmetric models consistently predict temperatures 1 to 2°C lower than the 
equivalent 3D results. As discussed in the benchmarking comparisons (Section 5.5), this is 
due to the increased recirculation associated with the axisymmetric jet. However, it is worth 
noting that the general trend and temperature profiles are consistent, and solutions can be 
obtained in under an hour. 


The circumferentially averaged 3D steady RANS model is in reasonable agreement with the 
URANS results, although some differences exist particularly around the centre of the CIP 
and along the IPT. This highlights the unsteady nature of the flow, although representative 
solutions can be generated over night. 


The circumferential and time-averaged URANS results provide the best match to the TALL- 
3D test data and WMLES model. The temperature profile on the IPT/ILW is consistent with 
the measurements, but the predictions are lower than expected. Although the IPT fluid only 
URANS predictions are lower than the WMLES results, the difference is not large enough 
to account for all of the differences in the IPT temperatures. This took 2 weeks to solve on 
192 compute cores, and so represents a significant increase in computational expense. This 
highlights the trade-off between accuracy and modelling approach that would be needed 
within any analysis process. 
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Validation Exercise (Full CHT) 


The baseline case temperatures are generally higher than the best estimate predictions and closer 
to the thermocouple measurements. This is due to the increased mass flow rate and heat loss ap- 
plied to the best estimate model. However, when the quantifiable values identified in the sensitivity 
analysis (Section 6.3.2) are calculated and compared to the TALL-3D test data (AT, = Tsim—Tmeas) 
in Table 6.7, the best estimate model significantly improves the outlet and external temperature 
predictions. 


Case Outlet External BP CIP IPT ILW Average 
Baseline case 1.7 23.4 0.7 0.2 -0.6 0.4 4.3 
Best estimate case 0.1 10.6 -0.4 -0.3 -2.0 -0.8 1.2 


Table 6.7: Summary of forced convection temperature difference results, AT;. 


It is worth noting that differences exist between the 3D URANS and WMLES results (Figure 6.8). 
However, although the WMLES results are slightly higher than the 3D URANS predictions, this 
difference is not enough to account for the differences between the full CHT predictions and the 
TALL-3D test data. 


This suggests that the best estimate model heat loss is more representative of the 3D test section 
conditions, but the heat loss may not be appropriately distributed in the 3D test section with the 
current best estimate parameters. This could be because the external heat loss from the lower part 
of the 3D test section is too high (i.e. below the CIP). If the heat loss through the top plate were 
increased even further, then the temperatures in the lower section would increase while maintaining 
the overall heat balance. This could be investigated by doing a more detailed automated calibration 
process on the mixed convection case. 


However, the low temperature predictions in the best estimate case may be due to an alternative 
reason, such as the relative location of the IPT and ILW thermocouples inside the 3D test section. 
As shown in Figure 6.9a, the IPT thermocouples are located inside a 6mm diameter tube that is 
open on one side, while the ILW thermocouples are located on the wall at the same azimuthal 
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Figure 6.9: Possible reason for differences between measurements and predictions. 


Under forced convection, the impinging jet will hit the bottom of the CIP and spread out radially. 
This is not relevant during mixed convection as the flow is stratified and the velocities in the 3D test 
section are small. This means that: 
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Validation Exercise (Full CHT) 


« Due to the position of the IPT tubes relative to the flow, a wake will be formed behind them. 
This slow moving fluid behind the IPT tubes is likely to be hotter than the surrounding areas 
due to the constant heat flux from the heater. 


« Although the four IPT tubes only represent 3% of the circumferential area, their location 
means that the wall temperature around the ILW thermocouples could be locally higher. 


An initial 3D RANS full CHT coarse model (10 million cells) with the IPT tubes and vanes included 
was solved steady-state and the results (Figure 6.9b) indicate that there could be a local heating 
of 1-2°C on the ILW behind the IPT tubes. However, in order to investigate this effect properly, 
the model would need to be solved unsteady with a more refined mesh. In addition, the simulation 
would have to solve for a much longer time period in order to ensure that the time-averaged tem- 
peratures at the four IPT locations (rather than circumferentially averaging the whole surface) are 
converged and stable. 


Another possible explanation could be the complexity of modelling the flow near the heated side 
wall. There, the decaying inertia of the radial flow reflected from the wall in the downward direction 
is fighting with buoyancy forces created by heating of the fluid as it moves along the wall. Note that 
the flow structure in the location of the IPTs is quite sensitive to the turbulence modelling approach 
(Figure 5.10). 


The differences between the baseline and best estimate case demonstrate the temperature varia- 
tion in the model due to the uncertainty in the boundary conditions, but also the need to consider 
and take into account simplifications that have been made in the modelling approach. In reality, the 
differences observed could be due to a combination of these effects. 


Despite some remaining unresolved differences, the fact that the fluid only 3D URANS predictions 
are in good agreement with the WMLES results, and the full CHT results are consistent with the 
TALL-3D measurements (within 2°C), demonstrates that the CFD model is capable of predicting 
forced convection in liquid metal and the importance of understanding and resolving any unsteady 
flow features. 
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7.1 


Application of the Results 


This case study has demonstrated the steps needed to apply and validate CFD modelling for 
forced and mixed convection in liquid metals with and without CHT. This has been achieved by 
undertaking a a simple hierarchical validation exercise using a heated channel basic test case 
(Section 2), followed by validation of forced and mixed convection predictions against TALL-3D 
measurement data (Section 4). 


This validation process has developed confidence in the ability of CFD tools to simulate forced and 
mixed convection in liquid metals, and provides evidence to apply them to similar thermal hydraulic 
phenomena as part of a design substantiation process. The lessons learnt and recommendations 
from this validation exercise are broken down into the following sections: 


* Section 7.1 provides general guidance for modelling liquid metal flows based on the heated 
channel and TALL-3D test cases. 


* The conclusions and recommendations from the TALL-3D validation exercise for mixed and 
forced convection flows are summarised in Sections 7.2 and 7.3 respectively. 


¢ Finally, Section 7.4 provides an indication of how this knowledge could be applied as part of 
a reactor design process. 


In addition, this case study has generated WMLES results using a simplified model of the TALL-3D 
test section under well-defined boundary conditions for mixed and forced convection scenarios. 
This provides benchmark simulations that allow a clear comparison and assessment of the accu- 
racy of RANS turbulence models. 


CFD modelling of Liquid Metal 


The heated channel and TALL-3D validation exercises have demonstrated that commercial CFD 
tools are capable of simulating forced and mixed convection flows in low Prandtl number fluids, but 
additional care needs to be taken when modelling liquid metal: 


* As with all CFD simulations, solution verification of the mesh and time step discretisation, 
convergence, etc. should be undertaken to ensure that the approach is appropriate for a 
particular application. 


* The sensitivity of the results of importance to different RANS turbulence models should be 
assessed, as this can vary depending on the specific application. Ideally, this should be 
undertaken by comparing the RANS simulations against experimental measurements or LES 
results of similar thermal hydraulic phenomena. 


¢ The sensitivity of the temperature predictions to modifying the turbulent Prandtl number 
should be investigated, as this can have a significant effect on the predicted temperatures. 
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7.2 


Application of the Results 


The heated channel study investigated the sensitivity of the temperature predictions to three mod- 
els for the turbulent Prandtl number (Section 2.2.2). These modifications represent adjustments to 
the turbulent heat transfer in existing RANS turbulence models for low Prandtl number fluids. As 
discussed in Volume 3 (Section 4.3) and Volume 5 (Section 4), more advanced turbulent heat flux 
models that are not dependent on the concept of turbulent Prandtl number, such as AHFM, are 
being developed. 


This is an active area of research that has the potential to significantly improve the accuracy of the 
RANS turbulence model simulations of low Prandtl number fluids. It is expected that over the next 
few years, these models will become more widely available within standard CFD analysis tools. 
Further validation will then be required to understand and assess their potential benefits for forced 
and mixed convection flows. 


Mixed Convection 


The mixed convection benchmark results (Section 5.4) suggested that for this test case a k- w SST 
turbulence model with full buoyancy terms and without the Kays modification provided the closest 
match to the WMLES results (within 1 to 2°C), as buoyancy effects are significant in this case. 


The low flow velocities resulted in thermal stratification in the 3D test section with an overturning 
negatively buoyant plume. The 3D CFD results show that large-scale unsteady effects are negligi- 
ble, and a 2D axisymmetric approach is appropriate. This significantly reduces the computational 
expense, and allows cases to be run in under an hour. 


Asensitivity analysis has been undertaken (Section 6.3.2) to investigate the impact of uncertainties 
in the material properties, inlet and boundary conditions on the predicted temperatures in the 3D 
test section. This confirmed that large variations in the temperature predictions can be achieved 
due to the sources of uncertainty in the experiment, and the uncertainty in the boundary conditions 
can be much larger than the uncertainty due to the modelling approach. For example, the baseline 
case predicted an outlet temperature that was almost 10°C higher than the measured value, as 
the external heat loss or heater power was inconsistent with the experimental test conditions. 


This finding emphasises the importance of conducting dedicated model calibration experiments 
where uncertainties in the boundary conditions can be properly quantified and possibly reduced. 


It is worth noting that, based on the mixed convection results, the mesh in the solid components 
could be significantly simplified and reduced by converting the air cavities to solid insulation, as 
their impact on the overall thermal resistance is small. In addition, the solid mesh could be coars- 
ened by using a non-conformal fluid-solid interface to reduce the solution times for industrial appli- 
cations. 


The mixed convection best estimate case demonstrates that a good level of agreement can be 
achieved with the model predictions within 2°C of all internal thermocouple measurements. 


¢ A significant reduction in the external thermal resistance was required to achieve the appro- 
priate heat loss. This demonstrates the large uncertainty in insulation thermal conductivity 
properties, and potential for under-estimation of the external HTC and/or redistribution of the 
heater power between the test section and outer insulation. 
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7.3 


Application of the Results 


* This demonstrates that an appropriate level of agreement has been achieved and that this 
modelling approach can be applied to simulate similar thermal hydraulic phenomena. 


For industrial applications, this analysis demonstrates the importance of understanding the uncer- 
tainty in the boundary conditions applied to CFD simulations, and assessing the impact of this 
uncertainty on the results of importance. Detailed information on the methods available to assess 
uncertainty is provided in Volume 4. 


Forced Convection 


The forced convection benchmark results (Section 5.5) demonstrated the complexity of modelling 
an impinging jet within a simple, cylindrical domain. Although the geometry is axisymmetric, large- 
scale unsteady fluctuations of the jet are predicted, which impact the level of recirculation within 
the 3D test section. 


For the forced convection case, a realizable k - « URANS turbulence model with the Kays modifi- 
cation provided the closest match to the WMLES results (< 1°C). The 2D axisymmetric results are 
consistent with the 3D URANS and WMLES predictions, but the temperatures are approximately 
1-2°C lower due to the increased recirculation below the CIP. 


* The 3D URANS results demonstrated that a good level of agreement has been achieved with 
the WMLES results and CFD is capable of simulating this type of forced convection flow. 


The full CHT forced convection results showed a similar trend between the 2D axisymmetric and 
3D URANS results, and highlighted the complexity of CHT simulations. In this case, the baseline 
case appears to be a better match to the internal thermocouples than the best estimate case, but 
the external heat loss is better predicted by the best estimate model. 


¢ This means that the heat loss may not be appropriately distributed in the 3D test section with 
the current best estimate parameters. This could be investigated and improved by doing a 
more detailed automated calibration process on the mixed convection case. 


Another possible reason for this difference could be because the detailed geometry of the 
IPT tubes has been neglected, which could cause a temperature increase local to the ILW 
thermocouples. This would require a significantly more detailed model to be developed and 
run for an extended period of time. 


In addition, the inherently unstable inertia driven downward flow in the vicinity of the heated 
wall that generates upward buoyancy forces is complex to model and could be sensitive to 
the turbulence modelling approach. 


From an industrial perspective, it is important to understand the trade-off between accuracy and 
turn-around i.e. 2D axisymmetric (one hour), 3D RANS (overnight) and 3D URANS (two weeks). 
Since the trends in the results between each method is consistent in this case, they could be used 
in combination to achieve a faster optimised design: 


* Once the relative accuracy of a 2D axisymmetric approach had been established, this method 
could be used to undertake a comparative design optimisation or uncertainty quantification 
approach. This would enable a large number of cases to be run quickly and enable the design 
space and parameters to be understood. 
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7.4 


Application of the Results 


* As the design solution matures, steady 3D RANS simulations could be used to confirm that 
the 2D axisymmetric results are representative and finalise the design parameters. 


¢ Finally, 3D URANS simulations of the final design and key uncertainty parameters could be 
undertaken to confirm the results of importance to an appropriate level of confidence. 


Reactor Application 


This case study forms part of a hierarchical validation process of forced and mixed convection in 
liquid metals, which could be used to support the safety justification for CFD analysis of reactor 
applications. The thermal hydraulic phenomena that have been validated by this analysis include: 


* Forced and mixed convection. 

* Thermal mixing. 

* Thermal stratification. 

* Propagation of a submerged jet. 
¢ Jet impingement on a surface. 


These thermal hydraulic phenomena are relevant to the upper and lower plenums of pool-type re- 
actors under normal operation (forced circulation) and loss of power (natural circulation) conditions. 
This case study could therefore be used to support the justification of the modelling approaches 
within a reactor application that involve these phenomena. However, it is important to ensure that 
the conditions within the reactor application are similar to the validation case by comparing the 
key non-dimensional parameters of the flow, such as Reynolds number, Péclet number, Rayleigh 
number and Richardson number. 


The validation results have demonstrated that CFD analysis can adequately resolve the fluid flow 
and heat transfer in liquid metal. In this case, the uncertainty in the boundary conditions has a 
more significant impact on the temperature predictions, compared to the differences in modelling 
approach. This highlights the importance of defining realistic boundary conditions and the benefit 
of uncertainty quantification to support the design substantiation within a graded approach. 


The next step in the validation process would be to apply these techniques to a more complex 
test facility, such as the E-SCAPE (European SCAled Pool Experiment) experimental facility. This 
thermal hydraulic test facility is a 1/6-scale model of the primary system of the MYRRHA (Multi- 
purpose hYbrid Research Reactor for High-tech Applications) reactor, with an electrical core sim- 
ulator, cooled by LBE. 


Alternatively, this method of generating higher fidelity LES predictions could be used to support 
the design of specific regions of a primary circuit, such as the lower plenum. Although it is not 
currently feasible to develop a LES model of a whole lower plenum, a smaller section with well- 
defined boundary conditions could be solved to provide confidence that a RANS model of the lower 
plenum flow could be used as part of the design process. 
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9 Nomenclature 


Latin Symbols 
A Area, m? 
At Atwood number (At = (p1 — p2)/(e1 + p2)) 
Bi Biot number (Bi = hL/kg) 
Cp, Cy Specific heat at constant pressure or volume, J kg~ K~* 
dorD_ Diameter (D, = 4A;;/Pcs for hydraulic diameter), m 
f Darcy-Weisbach friction factor 
Fo Fourier number (Fo = a,t/L?) 
Gr Grashof number (Gr = gL3Ap/v?p = gL?BAT/v?, using the Boussinesq approxi- 
mation Ap/p ~ —BAT, where AT is often taken as Ty, — Ts,00) 
g Acceleration due to gravity, ms~? 
h_ Specific enthalpy, J kg~+, Heat Transfer Coefficient (HTC), Wm~? K—! or height, m 
I Radiative intensity, W m~? sr~! or Wm~? sr~+ wm? for a spectral density, where sr 
(steradian) is solid angle 
J Radiosity, W m~? 
k Thermal conductivity, Wm-! K-? 
L Length or wall thickness, m 
M_ Molar mass of a species, kg kmol~? 
Ma _ Mach number (Ma = U/a, where a is the speed of sound) 
n_ Refractive index 
Nu Nusselt Number (Nu = AL/k¢) 
p Perimeter, m 
P Pressure (P; = static pressure, Py = total pressure), Nm~? or Pa 
Pe  Péclet number (Pe = RePr) 
Pr Prandtl number (Pr = cpu/ ke) 
q_ Heat flux (rate of heat transfer per unit area, g = Q/A), Wm~? 
Q_ Rate of heat transfer, W 
r Radius, m 
R_ Gas constant (for a particular gas, R = R/M), Jkg~!K~! 
R Universal gas constant, 8314.5 J kmol~! K~1 
Rip Thermal resistance, K W~! 
Ra Rayleigh number (Ra = GrPr) 
Re Reynolds number (Re = pUL/,, or for an internal flow Re = WD,,/A.<h) 
Ri Richardson number (Ri = Gr/Re’) 
Sr Strouhal number (Sr = fL/U, where f is frequency) 
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Nomenclature ‘Se, 


«x j<< 


Stanton number (St = Nu/RePr) 

Time, s 

Temperature (7, = static temperature, 77 = total temperature), K 
Wall friction velocity (u, = \/Tw/p), ms~+ 

Velocity, ms~? or thermal transmittance, Wm? K~! 
Specific volume, m* kg~? 

Volume, m? 

Mass flow rate, kg s~+ 
Wall distance, m 


Non-dimensional wall distance (y* = yu, /v) 


Greek Symbols 


Qn SS FE YM KF Hm R WR 


4 


oa 


Thermal diffusivity (a = k/pcp), m?s~? 


Volumetric thermal expansion coefficient (8 = —(1/)(0p/0T)), K~* 
Ratio of specific heats (y = cp/c,) 

Emissivity or surface roughness height, m 

Absorption coefficient, m~? 

Wavelength, m 


Viscosity, kg m~1s~? 


Kinematic viscosity and momentum diffusivity (v = w/p), m?s~? 
Density, kg m~? 

Stefan Boltzmann constant, 5.67 x 10-8 W m-? K~* 

Shear stress, Nm~? 


Porosity or void fraction 


Subscripts and Modifications 


b 
cs 
f 


i 


wn 


8 


~ 


Bulk (mass-averaged) quantity 
Cross-sectional quantity 

Quantity relating to a fluid 

Quantity relating to a particular species 
Total (stagnation) quantity 

Turbulent quantity 

Static quantity or quantity relating to a solid 
Quantity relating to a wall or surface 
Quantity far from a wall or in free-stream 
Average quantity 

Molar quantity 

Varying quantity 
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10 Abbreviations 


ADS 
AHFM 
Alsl 
BP 
CFD 
CFL 
CHT 
CIP 
DNS 
DO 
EVM 
FCH 
HTC 
IDDES 
ILW 
IPT 
KTH 
LBE 
LES 
NPP 
OLW 
R&D 
RANS 
RSM 
SA 
SBES 
SESAME 


SET 
SRQ 
SST 
TALL-3D 
THINS 
UDF 

UQ 
URANS 
VVUQ 
WMLES 


Accelerator-Driven System 

Algebraic Heat Flux Model 

American lron and Steel Institute 

Base Plate 

Computational Fluid Dynamics 

Courant-Friedrichs Lewy 

Conjugate Heat Transfer 

Circular Inner Plate 

Direct Numerical Simulation 

Discrete Ordinates 

Eddy Viscosity Model 

First Cell Height 

Heat Transfer Coefficient 

Improved Delayed Detached Eddy Simulation 
Inner Lateral Wall 

Inner Pipe Thermocouple 

Royal Institute of Technology 

Lead-Bismuth Eutectic 

Large Eddy Simulation 

Nuclear Power Plant 

Outer Lateral Wall 

Research and Development 

Reynolds-Averaged Navier-Stokes 

Reynolds Stress Model 

Sensitivity Analysis 

Stress-Blended Eddy Simulation 

Simulations and Experiments for the Safety Assessment of MEtal cooled reac- 
tors 

Separate Effect Test 

System Response Quantity 

Shear Stress Transport 

Thermal-hydraulic ADS Lead-bismuth Loop with 3D flow test section 
Thermal Hydraulics of Innovative Nuclear Systems 
User Defined Function 

Uncertainty Quantification 

Unsteady Reynolds-Averaged Navier-Stokes 
Verification, Validation and Uncertainty Quantification 
Wall Modeled Large Eddy Simulation 
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